Comparative transcriptome analysis of whiteflies raised on cotton leaf curl Multan virus-infected cotton plants

Cotton leaf curl Multan virus (CLCuMuV), a serious viral disease causative agent in cotton plants in South Asia, is transmitted by the Bemisia tabaci cryptic species complex in a persistent circulative manner. A previous study indicated that Asia II-7 whiteflies could transmit CLCuMuV, while Mediterranean (MED) whiteflies failed to transmit CLCuMuV. However, little is known about the genes involved in this process. In this study, Asia II-7 and MED B. tabaci were utilized to determine transcriptomic responses after 48 h of acquisition access periods (AAPs). Result of Illumina sequencing revealed that, 14,213 and 8,986 differentially expressed genes (DEGs) were identified. Furthermore, DEGs related to the immune system and metabolism of Asia II-7 and MED in response to CLCuMuV-infected plants were identified and analyzed using Gene Ontologies (GO) and the Kyoto Encyclopedia of Genes and Genomes (KEGG), and the number of related DEGs in MED was lower than that of Asia II-7. The most abundant groups of DEGs between both viruliferous and aviruliferous whitefly species were the zf-C2H2 family of transcription factors (TFs). Notably, in comparison to viruliferous MED, Asia II-7 exhibited more DEGs related to cathepsin biosynthesis. Overall, this study provides the basic information for investigating the molecular mechanism of how begomoviruses affect B. tabaci metabolism and immune response either as vector cryptic species or non-vector species.

1 Introduction

Host plant-mediated interaction between viruses and insects plays a vital role in the epidemiology of plant diseases because plant diseases caused by begomoviruses limit the production of many economic crops, such as tomato (1), cassava (2), and cotton (3), causing tremendous losses of agricultural crops worldwide (4, 5). Begomoviruses are transmitted in a persistent circulative manner by whiteflies Bemisia tabaci, where midgut and primary salivary glands have been identified as barriers in the circulation in the vector body (6–8). Overall, the interaction between plants, begomoviruses, and whiteflies is associated with various molecular mechanisms (9, 10).

Cotton leaf curl Multan virus (CLCuMuV) is a serious and economically important viral disease agent in cotton and ornamental plants in many regions of the world (11). In China, CLCuMuV was first detected in Guangzhou (12) and had rapidly spread to many provinces including Guangxi, Hainan, Fujian, and Jiangsu (13, 14). Among 39 known cryptic species of B. tabaci, two indigenous species, Asia II-1 and Asia II-7, are possibly the most efficient vectors of CLCuMuV (15, 16). On the contrary, two invasive species, Middle East Asia Minor 1 (MEAM1) and Mediterranean (MED), have been proven to be unable to transmit CLCuMuV. The previous study also revealed that after being ingested by whiteflies, CLCuMuV had a poorer capacity to cross the midguts of MEAM1 and MED than AsiaII1 (16). Additionally, when begomoviruses transmit from the midgut lumen to the hemolymph, they regulate numerous genes in whiteflies (17).

Previous reports compared the gene expression of whiteflies fed on healthy plants and those fed on virus-infected plants to unravel the mechanism of whiteflies and virus interaction (18–20). The expression of metabolism- and immune reaction-related genes of whiteflies was altered after the acquisition of Tomato yellow leaf curl China virus (TYLCCNV) and Tomato yellow leaf curl virus (TYLCV) (18, 21). However, to the best of our knowledge, few studies focused on the differences in gene expression between whitefly vector species and non-vector species, which could help reveal the interaction between begomovirus and whitefly.

In this study, we collected aviruliferous and viruliferous Asia II-7 and MED of B. tabaci, respectively, and compared all the genes to elucidate the interaction underlying different combinations of whiteflies and begomoviruses. Comparative results of aviruliferous/viruliferous Asia II-7 and aviruliferous/viruliferous MED were obtained. Overall, this study provides the theoretical foundation to explore the molecular mechanism of how begomoviruses affect the molecular response of whiteflies either as vectors or non-vectors and select candidate genes involved in the immune response and metabolism of whiteflies for future functional studies.

2 Materials and methods

2.1 Insects, plants, and viruses

Species of whiteflies, Asia II-7 and MED colony (mtCOI GeneBank accession code: KM821541 and GQ371165), and two cotton cultivars, namely, Gossypium hirsutum L. cv. Zhongmian40 and cv. 112–2, were used for experiments. Clones of CLCuMuV (accession number KP762786) and CLCuMuB (KP762787) were obtained from the plant protection institute, Guangdong Academy of Agricultural Sciences.

2.2 Preparation of whitefly samples

Cultures of whitefly Asia II-7 and MED were reared on uninfected cotton plants (G. hirsutum cv. Zhongmian40) in insect-proof cages separately. For experiments, uninfected cotton plants of Zhongmian40 were cultivated to 7–8 true leaf stages when used. For the preparation of CLCuMuV-infected plants, cotton plants of 112–2 were first cultivated to 2–3 true leaf stage when virus inoculation was conducted, and then, the virus-inoculated plants were further cultivated to 7–8 true leaf stage when used. The status of virus infection of these plants with typical symptoms was verified by PCR analysis, and the primers used for testing are presented in Supplementary Table S1. To obtain whiteflies fed on uninfected and CLCuMuV-infected plants, whitefly adults from Asia II-7 and MED cultures were collected 7–8 d post-emergence and then placed on uninfected and CLCuMuV-infected cotton plants of 112–2, respectively. The whiteflies were then collected after 48 h of AAP. For each of the four treatments, e.g., aviruliferous Asia II-7 (AA), viruliferous Asia II-7 (IU), aviruliferous MED (AM), and viruliferous MED (VM), three biological replicates in three separate cages were conducted. All whitefly cultures of the four treatments were reared in cages at 26 ± 1°C and 16 h light/8 h darkness.

2.3 Transcriptome sequencing and analysis

For each replicate sampling, 200 individual whiteflies were collected into 1.5 mL microcentrifuge tubes and then flash frozen using liquid nitrogen and stored at −80°C.

The total RNA of each sample was extracted using Trizol reagent (purchased from Takara, Japan), following the manufacturer’s instructions. Overall, 1.5 μg of RNA per sample was used as input material for RNA sample preparation. As described by Li et al. (17), RNA-Seq libraries were constructed using RNA extracted from whiteflies and sequenced by HiSeq 2,500 (Illumina, Inc. USA), and RNA-Seq data were analyzed as described by Ding et al. (22). RNA-Seq raw reads were processed and normalized to fragments per kilobase million mapped reads (FPKM), and differential expression analysis was performed using edger. The resulting p-values were adjusted for multiple testing using a false discovery rate (FDR). For the identification of differentially expressed genes (DEGs), the following cutoff parameters were used: the threshold of Q-value (adjusted p-value) of ≤0.001 and log₂ (fold change) ≥1. The deduced protein sequences of DEGs were used to generate Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) functional classifications. Data were analyzed by Novogene (Tianjin, China).

2.4 Transcriptome sequencing and analysis

Seven DEGs were selected for qPCR analysis to validate the data of RNA-Seq. Total RNA was extracted from Asia II-7 and MED whiteflies fed on either CLCuMuV-infected or uninfected cotton plants for 48 h using the same methods and kits of RNA-Seq. In total, 0.5 μg of total RNA was reversely transcribed using the HiScript 1st Strand cDNA Synthesis Kit (QIAGEN, USA). The qPCR was performed on the CFX96 TM Real-Time PCR Detection System (Bio-Rad, USA) with SYBR Premix ExTaq II (Takara, Japan). β-actin was used as the reference gene (20). The relative expression levels of DEGs were calculated by the 2^-ΔΔCt method (23). The qPCR primers of targeted and reference genes are presented in Table 1.

Table 1

Table 1. CLCuMuV-regulated candidate DEGs related to the CTS family in Asia II-7 and MED whiteflies.

2.5 CLCuMuV detection by polymerase chain reaction (PCR)

The DNA of the 10 whitefly adults of viruliferous and aviruliferous groups of Asia II-7 and MEAM1 was extracted with an Easy Pure Genomic DNA Kit (Trans Gen Biotech, Beijing, China), respectively. The DNA samples were stored at −20°C. The primer sequences used for CLCuMuV detection were forward (5’-CAGGAAGCAGGAAAATACGAGA-3′)/reverse (5’-TGGCAGTCCAACACAAAATACG-3′), with an amplicon of 831 bp. PCR was performed using 25 μL of a PCR mixture containing 2 μL of template DNA, 12.5 μL of 2xTaq Plus Master Mix (Vazyme, Nanjing, China), 1 μL of a 10 μM solution of each primer, and 8.5 μL of sterile double-distilled water. PCR amplification was performed with a GeneAmp 9,600 PCR system (Perkin–Elmer) using the following temperature program: 1 cycle of denaturation for 4 min at 94°C; 35 cycles of melting at 94°C for 45 s, annealing at 52°C for 45 s, and elongation at 72°C for 1 min; and a final extension round at 72°C for 10 min. The PCR amplicons were separated electrophoretically onto a 1% agarose gel (pulsed-field-certified agarose; Bio-Rad Laboratories AG).

3 Results

3.1 Transcriptome overview

The acquisition of CLCuMuV by Asia II-7 and MED was first confirmed by PCR (Fig. S1). Then, the transcriptome of cryptic species Asia II-7 and MED of B. tabaci fed on CLCuMuV-infected and healthy cotton plants was sequenced and compared. A total of 74.28 Gb of clean data and 117,542 unigenes were identified with annotation information, with a mean length of 952 bp of a gene. The N50 was 1,392 bp, and the Q30 (the probability of an incorrect base call 1 in 1,000 times) was more than 92.09% (Supplementary Table S2). A total of 39,977 genes showed a significant hit stnr_INV, 57.6% of which had a top hit to Bemisia tabaci followed by Nephila clavipes (4.8%) (Figure 1).

Figure 1

Figure 1. Nr homologous species distribution. All the genes were searched in the Nr nucleotide database using BLASTX.

3.2 DEGs in cryptic species Asia II-7 and MED of B. tabaci exposed to CLCuMuV

The number of DEGs in cryptic species Asia II-7 and MED of B. tabaci exposed to CLCuMuV infection was 14,213 and 8,986, respectively. The changed ratios of most DEGs in the two groups were between 2⁻⁵ and 2⁵ (Figure 2). Among all the 14,213 DEGs in viruliferous Asia II-7 and aviruliferous whiteflies, 8,735 were upregulated, while 5,478 were downregulated. As for the DEGs in viruliferous MED and aviruliferous whiteflies, 6,170 were upregulated, while 2,816 were downregulated. Moreover, aviruliferous/viruliferous Asia II-7 and aviruliferous/viruliferous MED shared 6,094 DEGs, with specific DEGs of 8,119 and 2,892, respectively (Figure 3).

Figure 2

Figure 2. Volcano plot of DEGs of AA vs. VA (A) and AM vs. VM (B). Non-DEGs are shown as blue dots. Upregulated and downregulated DEGs are shown as red and green dots, respectively. AA, aviruliferous Asia II-7; VA, viruliferous Asia II-7; AM, aviruliferous MED; and VM, viruliferous MED.

Figure 3

Figure 3. Venn diagram of DEGs of viruliferous/aviruliferous whiteflies. AA, aviruliferous Asia II-7; VA, viruliferous Asia II-7; AM, aviruliferous MED; and VM, viruliferous MED.

The GO analysis revealed that DEGs of aviruliferous/viruliferous Asia II-7 were associated with 34 GO items in three categories, i.e., biological process, molecular function, and cellular component (Figure 4). In the biological process category, dominant items were localization (865), the establishment of localization (843), transport (841), cellular process (864), and single-organism localization (679). In the molecular function category, transporter activity (517) and transmembrane transporter activity (453) were the dominant items. Moreover, in the cellular component category, the annotated genes were mostly involved in the membrane (1304) and membrane part (912). Finally, DEGs of aviruliferous/viruliferous MED were associated with 24 GO items, where dominant items were single-organism process (1434), membrane (922), membrane part (657), localization (624), the establishment of localization (612), and transport (608) (Figure 4).

Figure 4

Figure 4. Gene ontology analysis of AA vs. VA (A) and AM vs. VM (B). The bar chart shows the distribution of the corresponding GO terms (p < 0.05). Different colors represent different GO categories. AA, aviruliferous Asia II-7; VA, viruliferous Asia II-7; AM, aviruliferous MED; and VM, viruliferous MED.

In addition, KEGG analysis revealed that DEGs of aviruliferous/viruliferous Asia II-7 were mostly involved in lysosome (92), purine metabolism (68), antigen processing and presentation (58), protein digestion and absorption (54), starch and sucrose metabolism (52), and the cAMP signaling pathway (50) (Figure 5). The DEGs of aviruliferous/viruliferous MED were mostly enriched in the lysosome (69), purine metabolism (48), protein digestion and absorption (37), antigen processing and presentation (37), and the cAMP signaling pathway (37) (Figure 5).

Figure 5

Figure 5. KEGG pathway analysis of AA vs. VA (A) and AM vs. VM (B). Rich factor is the ratio of differentially expressed protein numbers annotated in this pathway. Greater rich factor (p < 1) means stronger intensiveness. p-value ranges from 0 to 1, and a lower p-value means greater intensiveness. Enriched pathway terms are displayed with p < 0.05. AA, aviruliferous Asia II-7; VA, viruliferous Asia II-7; AM, aviruliferous MED; and VM, viruliferous MED.

3.3 The highest enriched KEGG pathways of DEGs in VA vs. AA and VM vs. AM

The most abundant upregulated DEGs of VA vs. AA in the KEGG pathway were protein digestion and absorption (Supplementary Table S3), which involves a complex series of degradative processes by enzymes in the digestive system to break down ingested protein into amino acids and small peptides. These components are then transported into enterocytes via various transporters and then used for amino acid metabolism or protein synthesis. When comparing aviruliferous whiteflies with viruliferous ones, the enzymes marked in purple were upregulated in both Asia II-7 and MED after CLCuMuV infection, suggesting a common response to the virus in these two species. Enzymes marked in yellow showed upregulation only in Asia II-7 after infection (Figure 6A). Overall, collagen and KCN were only downregulated in AA vs. VA, indicating a species-specific response to CLCuMuV. Collagen is the most abundant protein in the body. Its fiber-like structure is used to make connective tissue and is a major component of the bone, skin, muscles, tendons, and cartilage. The significant upregulation of collagen in viruliferous Asia II-7 indicated that the CLCuMuV might be impacting the connective tissue structure of the host and repair mechanisms of Asia II-7. “KCN” refers to a family of genes encoding potassium channels, which are crucial for maintaining the electrochemical gradient across cell membranes. Thus, upregulation of the KCN process in VA vs. AA might be a compensatory response to maintain cellular function, which was disturbed by CLCuMuV.

Figure 6

Figure 6. The highest enriched KEGG pathways of DEGs in VA vs. AA. (A) Upregulation DEGs; (B) downregulation DEGs. AA, aviruliferous Asia II-7; and VA, viruliferous Asia II-7. The purple or green box represents upregulation or downregulation in both viruliferous Asia II-7 and MED. The yellow or red box represents upregulation or downregulation only in viruliferous Asia II-7.

The most abundant downregulated DEGs of VA vs. AA in the KEGG pathway were ribosomes, which show the composition of ribosomes (Figure 6B). To our surprise, 44 ribosome compositions were downregulated in viruliferous Asia II, with only two of them (L37Ae and S21e) shared by both viruliferous Asia II and MED. The CLCuMuV impacted practically all ribosome biosynthesis activities, including initiation factors (IFs), elongation factors (EFs), release factors (RFs), and RNA polymerase subunits (Rpo), showing a deleterious impact on the ribosomes of Asia II.

3.4 Specific DEGs in Asia II-7and MED whitefly fed on CLCuMuV-infected cotton plants

A total of 8,119 specific DEGs in Asia II-7 MED whiteflies fed on CLCuMuV-infected cotton plants were identified, among which 4,342 genes were upregulated and 3,777 were downregulated. A total of 30 GO terms with the most significant enrichment were selected (Figure 7). In the biological process category, eight and nine items were assigned to upregulated and downregulated DEGs, respectively. Specifically, proteolysis, vitamin metabolic process, and hexachlorocyclohexane metabolic process were assigned to upregulated DEGs, while peptide metabolic, translation, and macromolecule biosynthetic processes were assigned to downregulated DEGs. As for the molecular function category, 22 items including transferring hexosyl groups, cysteine-type peptidase activity, and transferring glycosyl groups were associated with upregulated DEGs. On the contrary, seven items such as structural constituent of activity, structural molecule activity, and metal ion binding and transporter activity were assigned to downregulated DEGs. Finally, in the cellular component category, 12 items including plasma membrane, cell, and ribosome were linked with downregulated DEGs.

Figure 7

Figure 7. Specific DEGs in viruliferous Asia II-7 whiteflies (A, upregulated; B, downregulated) and in viruliferous Asia II-7 whiteflies (C, upregulated; D, downregulated) by gene ontology analysis. The bar chart shows the distribution of corresponding GO terms (p < 0.05). Different colors represent different GO categories.

For MED whiteflies fed on CLCuMuV-infected cotton plants, a total of 2,892 specific DEGs were identified with 1,778 upregulated genes and 1,114 downregulated genes. In total, 30 GO terms with the most significant enrichment were selected (Figure 7). In the biological process category, 12 and 20 items were assigned to upregulated and downregulated DEGs, respectively. Oxidation–reduction process, metabolic process, and electron transport chain were assigned to upregulated DEGs, while sodium ion export across the plasma membrane, proline metabolic process, and export from cells were assigned to downregulated DEGs. In the molecular function category, upregulated DEGs were found in 13 items, such as oxidoreductase activity, cofactor binding, and hemebinding. On the contrary, downregulated DEGs were associated with 10 items, such as oxaloacetate decarboxylase activity, transporter activity, and carboxylase activity. As for the cellular component category, three items including respiratory chain complex, cytochrome complex, and respirasome were assigned to upregulated DEGs, whereas downregulated DEGs were linked with membrane.

3.5 Common response of Asia II-7 and MED whiteflies fed on CLCuMuV-infected cotton plants

A total of 6,094 common DEGs in Asia II-7 and MED whiteflies fed on CLCuMuV-infected plants were found. A total of 30 GO terms with the most significant enrichment were selected. As a result, most of the common DEGs in whiteflies either as vector species or non-vector species were associated with the membrane (Figure 8). In the biological process category, upregulated DEGs were found in 10 items including carbohydrate metabolic process, biological process, and transmembrane transport. Downregulated DEGs were annotated to four items, including cellular response to oxidative process and transport process, amide transport, and cation transport. In the molecular function category, 15 items including transporter activity, G protein-coupled receptor activity, and molecular function were assigned to upregulated DEGs. For comparison, downregulated DEGs were associated with five processes, including aspartic-type endopeptidase activity, aspartic-type peptidase activity, and cation transmembrane transporter activity. In the cellular component category, upregulated DEGs were assigned to 4 items, such as an intrinsic component of membrane, an integral component of membrane, and membrane, whereas downregulated DEGs were assigned to 12 items including the membrane part, cellular components, and the organelle part.

Figure 8

Figure 8. Common Response of Asia II-7 and MED whiteflies fed on CLCuMuV-infected cotton plants by gene ontology analysis. The bar chart shows the distribution of corresponding GO terms (p < 0.05). Different colors represent different GO categories. Left, upregulated; Right, downregulated.

3.6 DEGs between VA vs. AA and VM vs. AM in cathepsin (CTS) biosynthesis

The CTS family of proteases plays a crucial role in the immune response, particularly through the process of antigen processing and presentation. We compared the DEGs between VA vs. AA and VM vs. AM in the CTS family (Table 1). Overall, 45 DEGs in the CTS family were identified in VA vs. AA, with 42 genes upregulated in viruliferous Asia II-7. The number of DEGs in the CTS family in VM vs. AM was much less, with 18 upregulated genes and 2 downregulated genes. Most DEGs were identified to belong to the CTSB family in both whitefly species.

3.7 The most abundant transcription factors between Asia II-7 and MED whiteflies induced by CLCuMuV infection

We further investigated and compared the TFs-related variations between Asia II-7 and MED whiteflies after CLCuMuV infection. Overall, 204 and 121 differentially expressed TFs, divided into 37 and 28 families, were observed in VA vs. AA and VM vs. AM, respectively (Supplementary Table S4). Most TFs in both whiteflies exhibited upregulated after CLCuMuV infection, with a proportion of 77.9% II-7 and 87.6% MED in Asia. The families with more than three DEGs are presented in Table 2. The zf-C2H2 TF family showed the largest number of DEGs, with 70 genes differentially expressed; of which, 37 were upregulated and 33 were downregulated in the AA vs. VA, and 34 DEGs were identified in the AM vs. VM. The TF families including homeobox, ZBTB, Forkhead, bZIP, and bHLH exhibited notable response patterns to virus infection, with more than 10 DEGs in AA vs. VA.

Table 2

Table 2. The most abundantly detected differentially expressed TFs among whiteflies in response to CLCuMuV infection.

3.8 RT-qPCR validation of DEGS

The differential expression of seven selected DEGs was verified by RT-qPCR. As a result, five genes were upregulated in Asia II-7 whiteflies fed on CLCuMuV-infected cotton plants, while two genes were downregulated in Asia II-7 whiteflies. All these candidate genes were related to the immune and single-organism processes in aviruliferous/viruliferous Asia II-7 and aviruliferous/viruliferous MED. Specifically, Cluster-35353.31682 (RHOBTB1), Cluster-35353.34860 (SCD), Cluster-35353.19909 (CTSB), and Cluster-35353.38911(ACTB_G1) were important for the transport of B. tabaci, the activity of metabolism, and immune pathways. Other genes such as Cluster-35353.38083 (HK), Cluster-35353.35641(HDAC4_5), and Cluster-11110.0 (EIF4E) were associated with translation activity. Moreover, all seven genes showed the same expression pattern (upregulated or downregulated) for both RT-qPCR and RNA-Seq analysis (Figure 9; Supplementary Table S5).

Figure 9

Figure 9. Expression patterns generated by the results of RT-qPCR. Five upregulated and two downregulated DEGs identified by RNA-Seq were validated by comparative CT (2^-ΔΔCT) RT-qPCR with β-actin as the reference gene. All data were subjected to one-way analysis of variance (ANOVA) to determine the significant differences between aviruliferous and viruliferous whiteflies. Different small letters are significantly different at a probability level of 0.05.

4 Discussion

To investigate the interaction between begomoviruses and vector and non-vector cryptic species of whiteflies, gene expression profiles of viruliferous Asia II-7 (can transmit CLCuMuV) and MED (cannot transmit CLCuMuV) whiteflies were demonstrated in this study. Our results suggested the number of DEGs in Asia II-7 exposed to CLCuMuV (14,213) was higher than that exposed to MED (8,916) (Figure 4). The increased number of DEGs might indicate that Asia II-7 whiteflies activate a broader array of genes in response to the virus, to manage the infection or mediate transmission processes. Similarly, a study by Li et al. (17) suggested that the number of DEGs in Sitobion avenae exposed to barley yellow dwarf virus (BYDV) (a virus vectored by S. avenae) was higher than that exposed to wheat dwarf virus (WDV) (a virus not vectored by S. avenae). These results indicated that plant viruses might indirectly affect vectors compared with the effect on non-vectors. However, it still remains uncertain whether the increased number of DEGs was due to a stronger impact of CLCuMuV on Asia II-7 whiteflies.

Researchers have demonstrated that begomovirus infection induced transcriptional changes in the genes in whiteflies (17, 20, 21). Previous reports found that genes clustered in lipid, carbohydrate, and amino acid metabolism were downregulated in TYLCV-infected whiteflies (17, 18). It was also suggested that plant viruses could induce the fecundity and longevity of whiteflies, which might alter the expression patterns of genes involved in metabolic processes (19, 24, 25). In this study, GO enrichment analysis showed that the DEGs in viruliferous Asia II-7 and MED whiteflies were all highly enriched in metabolic process, transport, translation, and cytoskeleton (Figures 7, 8). This suggests that begomovirus infection influences similar metabolic pathways in different whitefly populations, supporting previous findings while also highlighting specific pathways that impacted our study.

Immune-related genes were important for resistance to begomovirus infection in whiteflies (26). Immune systems also impact virus replication and transmission by vectors. The previous study revealed that several immune pathways are vital for insects against viral infections, including Janus kinase, signal transducer and activator of transcription (JAK/STAT), Toll (27), immune deficiency, and c-Jun N-terminal kinase (28), and RNAi (29). In this study, when compared with MED whiteflies fed on viruliferous cotton plants, more immune-related genes were found upregulated in CLCuMuV-Asia II-7, such as cathepsin (CTS) family (Figure 6). Specifically, our data revealed that 45 DEGs in the CTS family were significantly expressed in VA vs. AA, with 42 of these genes upregulated. This contrasts with only 20 DEGs identified in the CTS family in VM vs. AM, highlighting a more robust immune response in Asia II-7 whiteflies. The increased expression of CTSB in Asia II-7 compared with MED was further validated by qPCR analysis (Figure 9), emphasizing the significant immune activation in this strain. The upregulation of CTS genes in our study indicates a potential mechanism through which CLCuMuV influences the immune system of whitefly. Li et al. (30) identified many DEGs associated with immune response, with five genes (including CTSB) upregulated in S. avenae-exposed BYDV. This alignment with our findings highlights the relevance of CTS genes in the immune response to viral infections. In viral infections, CTS can facilitate virus entry, replication, and spread by activating viral glycoproteins, processing antigens, and triggering cell death (31). There have been studies indicating an increase in the CTS families in whiteflies, which is linked to the whitefly’s immunological responses during the process of acquiring viruses or other interactions between whitefly and viruses (32). The dramatic upregulation of CTS genes in viruliferous Asia II-7 suggests that CLCuMuV has a significant detrimental impact on this whitefly population, potentially through heightened immune activation leading to increased metabolic costs or other deleterious effects. Viral transportation plays an important role in virus transmission, and whitefly cellular proteins affect the interaction between whiteflies and viruses (17). Whiteflies are non-vectors of some begomoviruses because of the presence of a selective transmission barrier at the luminal membrane surface of epithelial cells of the midgut (16, 33). Cyclophilin B and collagen (34), sugar transporter ERD6-like 6, ATP-binding cassette (ABC) transporters, and other proteins have been reported to influence viral transport (18, 20, 35). In this study, compared with AM vs. VM, specific upregulated DEGs in AA vs. VA were assigned to catalytic activity, binding, and transporter activity. Specifically, the expression of genes related to potassium channel protein was significantly upregulated in CLCuMuV- Asia II-7 (Figure 6). The subfamilies in KCN, such as KCNQ, KCNH, and KCNN, are important for the nervous system, cardiac tissue, and T-cell activation. It has been reported that the replication of the Israeli acute paralysis virus (IAPV) in bees can be regulated by potassium channel expression (36). Our data revealed that the expression of several potassium channel proteins was significantly elevated in CLCuMuV-infected Asia II-7 compared with their non-viruliferous counterparts. This upregulation suggests an adaptive response, where increased potassium channel activity could help mitigate cellular stress caused by viral infection. These proteins might influence the ionic balance and cell signaling pathways, which are critical during viral infection and replication. By contrast, the expression of ribosome components decreased significantly in VA vs. AA (Figure 6), which could be due to a strategy by the virus to suppress the host’s protein synthesis machinery, thereby reducing the ability of the host to mount an effective immune response (37). The significant downregulation of ribosomal proteins in viruliferous Asia II-7 whiteflies suggests that CLCuMuV might be manipulating the host’s translational machinery to favor its replication and persistence. By reducing the capacity of the host to produce proteins, the virus might be evading immune detection and response. Interestingly, the KCN and ribosome were not significantly affected in CLCuMuV-MED. These findings may help explain why MED has a lower viral transmission rate than Asia II-7.

Transcription factors are crucial regulatory proteins that control the expression of genes, including those involved in immune responses, cell cycle regulation, and stress responses. Moreover, our study focused on the TF regulated by CLCuMuV infection. Most differentially expressed TFs were upregulated in viruliferous whiteflies, with C2H2 ZF being the most abundant TFs (Table 2). The results were consistent with the previous report that the most abundant TFs among viruliferous whiteflies belonged to the zf-C2H2 family (38). It has been reported that ZAD containing C2H2 ZF proteins (ZAD-ZNF) constitute the most abundant class of insect transcription factors (39). Silencing of the ZF-like gene significantly increases virus loads retained within the whiteflies using three diverse begomoviruses (40). In our study, the upregulation of C2H2 ZF transcription factors in CLCuMuV-infected whiteflies suggests a critical role in the response of whiteflies to viral infection. Specifically, we identified that the expressions of multiple C2H2 ZF genes were markedly higher in viruliferous Asia II-7 whiteflies compared to their non-viruliferous counterparts. These transcription factors are likely responsible for controlling a group of genes that improve the defensive mechanisms of whitefly against the virus. This might possibly restrict the reproduction and dissemination of the virus inside the host. The substantial increase in the expression of these transcription factors may suggest a strong transcriptional response aimed at stimulating immunological pathways and other defensive mechanisms, enhancing resistance to the virus. Homeobox TFs were also abundant in whiteflies (Table 2). These genes are recognized for encoding proteins that include a homeodomain, which has significant functions in several developmental processes such as segmentation, differentiation, and control of the cell cycle. In fruit flies, the homeodomain protein POU regulates the continuous expression of antimicrobial peptide (AMP) genes, together with other regulators, to enhance the first defense against infection (41, 42). In this investigation, we have discovered 23 homeobox transcription factors that are expressed differently between VA and AA (Table 2). This finding suggests that homeobox genes play a role in the response to viruses. The differential expression of these homeobox genes suggests that they may be involved in regulating genes critical for antiviral responses, possibly through pathways that control immune function and stress responses. In the present study, we also identified 10 and 9 differently expressed bZIPs in VA vs. AA and VM vs. AM, respectively. BZIPs have significant functions in several biological processes in plants, such as development phases and the response to both abiotic and biotic stressors and signals (43). Nevertheless, the precise function of the bZIPs in insects remains unclear. Therefore, it is essential to conduct a functional analysis of whitefly transcription factors and the genes they control in future studies in order to get a comprehensive knowledge of the cellular processes involved in viral transmission.

5 Conclusion

In summary, we uncovered transcriptional responses of Asia II-7 and MED B. tabaci after feeding CLCuMuV-infected plants for 48 h. Our results indicated that CLCuMuV might impact vector cryptic species more than non-vector cryptic species because CLCuMuV-Asia II-7 had more DEGs (14,213) than CLCuMuV-MED (8,916). Genes involved in metabolism and immune systems were upregulated in viruliferous whiteflies. The most abundant upregulated DEGs of VA vs. AA in the KEGG pathway were protein digestion and absorption, while the most enriched downregulated pathway was ribosome. Furthermore, we identified specific genes in whiteflies related to cathepsin biosynthesis and transcription factors after virus infection. Overall, this study lays the foundation for further study on whiteflies and begomovirus interactions.

Data availability statement

The data presented in the study are deposited in the NCBI repository, the SRA accession number is PRJNA1144620.

Author contributions

TC: Writing – original draft, Writing – review & editing, Data curation, Formal analysis, Methodology. YJ: Data curation, Writing – review & editing. JC: Methodology, Writing – review & editing. GQ: Formal analysis, Investigation, Project administration, Writing – review & editing.

Funding

The author(s) declare that financial support was received for the research, authorship, and/or publication of this article. This study was funded by the National Natural Science Foundation of China (32001973), the National Key R&D Program of China (2022YFC2601400), and the Innovation Team of Modern Agricultural Industry Technology System of Guangdong Province (2023KJ122).

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.

Supplementary material

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

References

1. Scholthof, KBG, Adkins, S, Czosnek, H, Palukaitis, P, Jacquot, E, Hohn, T, et al. Top 10 plant viruses in molecular plant pathology. Mol Plant Pathol. (2011) 12:938–54. doi: 10.1111/j.1364-3703.2011.00752.x

PubMed Abstract | Crossref Full Text | Google Scholar

2. Legg, JP, Shirima, R, Tajebe, LS, Guastella, D, Boniface, S, Jeremiah, S, et al. Biology and management of Bemisia whitefly vectors of cassava virus pandemics in Africa. Pest Manag Sci. (2014) 70:1446–53. doi: 10.1002/ps.3793

PubMed Abstract | Crossref Full Text | Google Scholar

3. Navas-Castillo, J, Fiallo-Olive, E, and Sanchez-Campos, S. Emerging virus diseases transmitted by whiteflies. Annu Rev Phytopathol. (2011) 49:219–48. doi: 10.1146/annurev-phyto-072910-095235

Crossref Full Text | Google Scholar

4. Hogenhout, SA, Ammar, E-D, Desouky, A, Whitfield, AE, and Redinbaugh, MG. Insect vector interactions with persistently transmitted viruses. Annu Rev Phytopathol. (2008) 46:327–59. doi: 10.1146/annurev.phyto.022508.092135

Crossref Full Text | Google Scholar

5. Ye, J, Gong, YQ, and Fang, RX. Research progress and perspective of tripartite interaction of virus-vector-plant in vector-borne viral diseases. Bull Chin Acad Sci. (2017) 32:845–55. doi: 10.16418/j.issn.1000-3045.2017.08.006

Crossref Full Text | Google Scholar

6. Guo, T, Zhao, J, Pan, LL, Geng, L, Lei, T, Wang, XW, et al. The level of midgut penetration of two begomoviruses affects their acquisition and transmission by two species of Bemisia tabaci. Virology. (2018) 515:66–73. doi: 10.1016/j.virol.2017.12.004

PubMed Abstract | Crossref Full Text | Google Scholar

7. Pan, LL, Chen, QF, Zhao, JJ, Guo, T, Wang, XW, Hariton-Shalev, A, et al. Clathrin-mediated endocytosis is involved in tomato yellow leaf curlvirus transport across the midgut barrier of its whitefly vector. Virology. (2017) 502:152–9. doi: 10.1016/j.virol.2016.12.029

PubMed Abstract | Crossref Full Text | Google Scholar

8. Wei, J, Zhao, JJ, Zhang, T, Li, FF, Ghanim, M, Zhou, XP, et al. Specific cells in the primary salivary glands of the whitefly Bemisia tabaci control retention and transmission of begomoviruses. J Virol. (2014) 88:13460–8. doi: 10.1128/jvi.02179-14

PubMed Abstract | Crossref Full Text | Google Scholar

9. Götz, M, Popovski, S, Kollenberg, M, Gorovits, R, Brown, JK, Cicero, JM, et al. Implication of Bemisiatabaci heat shock protein 70 in begomovirus-whitefly interactions. J Virol. (2012) 86:13241–52. doi: 10.1128/JVI.00880-12

PubMed Abstract | Crossref Full Text | Google Scholar

10. Rana, VS, Popli, S, Saurav, GK, Raina, HS, Chaubey, R, Ramamurthy, VV, et al. A Bemisia tabaci midgut protein interacts with begomoviruses and plays a role in virus transmission. Cell Microbiol. (2015) 18:663–78. doi: 10.1111/cmi.12538

Crossref Full Text | Google Scholar

11. Mansoor, S, Briddon, RW, Bull, SE, Bedford, ID, Bashir, A, Hussain, M, et al. Cotton leaf curl disease is associated with multiple monopartite begomoviruses supported by a single DNA β. Arch Virol. (2003) 148:1969–86. doi: 10.1007/s00705-003-0149-y

PubMed Abstract | Crossref Full Text | Google Scholar

12. Mao, MJ, He, ZH, Yu, H, Chen, BS, Meng, JR, and Qin, BX. Molecular characterization of cotton leaf curl Multan virus and its satellite DNA that infects Hibiscus rosa-sinensis. Chin J Virol. (2008) 24:64–8. doi: 10.13242/j.cnkibingduxuebao.001881

Crossref Full Text | Google Scholar

13. Cai, JH, Xie, K, Lin, L, Qin, BX, Chen, BS, Meng, JR, et al. Cotton leaf curl Multan virus newly reported to be associated with cotton leaf curl disease in China. J Plant Pathol. (2010) 59:794–5. doi: 10.1111/j.1365-3059.2010.02266.x

Crossref Full Text | Google Scholar

14. Du, ZG, Tang, YF, He, ZF, and She, XM. High genetic homogeneity points to a single introduction event responsible for invasion of cotton leaf curl Multan virus and its associated betasatellite into China. Virology. (2015) 12:163. doi: 10.1186/s12985-015-0397-y

PubMed Abstract | Crossref Full Text | Google Scholar

15. Chen, T, Saeed, Q, He, ZF, and Lu, L. Transmission efficiency of cotton leaf curl Multan virus by three cryptic species of Bemisia tabaci complex in cotton cultivars. PeerJ. (2019) 7:e7788. doi: 10.7717/peerj.7788

PubMed Abstract | Crossref Full Text | Google Scholar

16. Pan, LL, Cui, XY, Chen, QF, Wang, XW, and Liu, SS. Cotton leaf curl disease: which whitefly is the vector? Phytopathology. (2018) 108:1172–83. doi: 10.1094/PHYTO-01-18-0015-R

PubMed Abstract | Crossref Full Text | Google Scholar

17. Li, M, Zhao, J, and Su, YL. Transcriptome analysis of gene expression profiles of tomato yellow leaf curl virus-infected whiteflies over different viral acquisition access periods. Insects. (2020) 11:297. doi: 10.3390/insects11050297

PubMed Abstract | Crossref Full Text | Google Scholar

18. Hasegawa, DK, Chen, W, Zheng, Y, Kaur, N, Wintermantel, WM, Simmons, AM, et al. Comparative transcriptome analysis reveals networks of genes activated in the whitefly, Bemisia tabaci, when fed on tomato plants infected with tomato yellow leaf curl virus. Virology. (2018) 513:52–64. doi: 10.1016/j.virol.2017.10.008

PubMed Abstract | Crossref Full Text | Google Scholar

19. Luan, JB, Li, JM, Varela, N, Wang, YL, Li, FF, Bao, YY, et al. Global analysis of the transcriptional response of whitefly to tomato yellow leaf curl China virus reveals the relationship of coevolved adaptations. J Virol. (2011) 85:3330–40. doi: 10.1128/JVI.02507-10

PubMed Abstract | Crossref Full Text | Google Scholar

20. Zhao, J, Chi, Y, Zhang, XJ, Lei, T, Wang, XW, and Liu, SS. Comparative proteomic analysis provides new insight into differential transmission of two begomoviruses by a whitefly. Virology. (2019) 16:32. doi: 10.1186/s12985-019-1138-4

PubMed Abstract | Crossref Full Text | Google Scholar

21. Li, M, Hu, J, Xu, F, and Liu, SS. Transmission of Tomato yellow leaf curl virus by two invasive biotypes and a Chinese indigenous biotype of the whitefly Bemisia tabaci. Int J Pest Manage. (2010) 56:275–80. doi: 10.1080/09670871003743428

PubMed Abstract | Crossref Full Text | Google Scholar

22. Ding, TB, Chen, EH, Niu, JZ, and Chu, D. Transcriptome profiling of the whitefly Bemisia tabaci MED in response to single infection of tomato yellow leaf curl virus, tomato chlorosis virus, and their co-infection. Front Physiol. (2019) 10:302. doi: 10.3389/fphys.2019.00302

PubMed Abstract | Crossref Full Text | Google Scholar

23. Livak, KJ, and Schmittgen, TD. Analysis of relative gene expression data using real time quantitative PCR and the 2^(–ΔΔCT)method. Methods. (2001) 25:402–8. doi: 10.1006/meth.2001.1262

Crossref Full Text | Google Scholar

24. Pusag, JCA, Hemayet Jahan, SM, Lee, KS, Lee, S, and Lee, KY. Upregulation of temperature susceptibility in Bemisia tabaci upon acquisition of tomato yellow leaf curl virus (TYLCV). J Insect Physiol. (2012) 58:1343–8. doi: 10.1016/j.jinsphys.2012.07.008

PubMed Abstract | Crossref Full Text | Google Scholar

25. Rubinstein, G, and Czosnek, H. Long-term association of tomato yellow leaf curl virus with its whitefly vector Bemisia tabaci: effect on the insect transmission capacity, longevity and fecundity. J Gen Virol. (1997) 78:2683–9. doi: 10.1099/0022-1317-78-10-2683

PubMed Abstract | Crossref Full Text | Google Scholar

26. Wang, ZZ, Shi, M, Huang, YC, Wang, XW, Stanley, D, and Chen, XX. A peptidoglycan recognition protein acts in whitefly (Bemisia tabaci) immunity and involves in begomovirus acquisition. Sci Rep. (2016) 6:37806. doi: 10.1038/srep3780

PubMed Abstract | Crossref Full Text | Google Scholar

27. Lemaitre, B, and Hoffmann, J. The host defense of Drosophila melanogaster. Annu Rev Immunol. (2007) 25:697–743. doi: 10.1146/annurev.immunol.25.022106.141615

Crossref Full Text | Google Scholar

28. Zhang, Z, Zhang, P, Li, W, Zhang, J, Huang, F, Yang, J, et al. De novo transcriptome sequencing in Frankliniella occidentalis to identify genes involved in plant virus transmission and insecticide resistance. Genomics. (2013) 101:296–305. doi: 10.1016/j.ygeno.2013.02.005

PubMed Abstract | Crossref Full Text | Google Scholar

29. Wang, XH, Aliyari, R, Li, WX, Li, HW, Kim, K, Carthew, R, et al. RNA interference directs innate immunity against viruses in adult drosophila. Science. (2006) 312:452–4. doi: 10.1126/science.1125694

PubMed Abstract | Crossref Full Text | Google Scholar

30. Li, DD, Su, D, Tong, ZQ, Zhang, C, Zhang, G, Zhao, H, et al. Virus-dependent and-independent responses of Sitobion avenae (Homoptera: Aphididae) feeding on wheat infected by transmitted and nontransmitted viruses at transcriptomic level. J Econ Entomol. (2019) 112, 2067–2076. doi: 10.1093/jee/toz162

PubMed Abstract | Crossref Full Text | Google Scholar

31. Scarcella, M, d’Angelo, D, Ciampa, M, Tafuri, S, Avallone, L, Pavone, LM, et al. The key role of lysosomal protease cathepsins in viral infections. Int J Mol Sci. (2022) 23:9089. doi: 10.3390/ijms23169089

PubMed Abstract | Crossref Full Text | Google Scholar

32. Chen, W, Hasegawa, DK, Kaur, N, Kliot, A, Pinheiro, PV, Luan, J, et al. The draft genome of whitefly Bemisia tabaci MEAM1, aglobal crop pest, provides novel insights into virus transmission, hostadaptation, and insecticide resistance. BMC Biol. (2016) 14:110. doi: 10.1186/s12915-016-0321-y

PubMed Abstract | Crossref Full Text | Google Scholar

33. Ohnishi, J, Kitamura, T, Terami, F, and Honda, K. A selective barrier in the midgut epithelial cell membrane of the non-vector whitefly Trialeurodes Vaporariorum to tomato yellow leaf curl virus uptake. J Gen Plant Pahtol. (2009) 75:131–9. doi: 10.1007/s10327-009-0147-3

Crossref Full Text | Google Scholar

34. Kanakala, S, and Ghanim, M. Implication of the whitefly Bemisia tabaci cyclophilin B protein in the transmission of tomato yellow leaf curl virus. Front Plant Sci. (2016) 7:1702. doi: 10.3389/fpls.2016.01702

PubMed Abstract | Crossref Full Text | Google Scholar

35. Geng, L, Qian, LX, Shao, RX, Liu, YQ, Liu, SS, and Wang, XW. Transcriptome profiling of whitefly guts inresponse to tomato yellow leaf curl virus infection. Virol J. (2018) 15:14. doi: 10.1186/s12985-018-0926-6

PubMed Abstract | Crossref Full Text | Google Scholar

36. Fellows, CJ, Simone-Finstrom, M, Anderson, TD, and Swale, DR. Potassium ion channels as a molecular target to reduce virus infection and mortality of honey bee colonies. Virol J. (2023) 20:134. doi: 10.1186/s12985-023-02104-0

PubMed Abstract | Crossref Full Text | Google Scholar

37. Chen, S, Fu, Z, Chen, K, Zheng, X, and Fu, Z. Decoding HiPSC-CM’s response to SARS-CoV-2: mapping the molecular landscape of cardiac injury. BMC Genomics. (2024) 25:271. doi: 10.1186/s12864-024-10194-5

PubMed Abstract | Crossref Full Text | Google Scholar

38. Farooq, T, Lin, Q, She, X, Chen, T, Tang, Y, and He, Z. Comparative transcriptome profiling reveals a network of differentially expressed genes in Asia II 7 and MEAM1 whitefly cryptic species in response to early infection of cotton leaf curl Multan virus. Front Microbiol. (2022) 13:1004513. doi: 10.3389/fmicb.2022.1004513

PubMed Abstract | Crossref Full Text | Google Scholar

39. Krystel, J, and Ayyanathan, K. Global analysis of target genes of 21 members of the ZAD transcription factor family in Drosophila melanogaster. Gene. (2013) 512:373–82. doi: 10.1016/j.gene.2012.09.114

PubMed Abstract | Crossref Full Text | Google Scholar

40. Ghosh, S, Srinivasan, R, and Ghanim, M. A C2H2 zinc finger transcription factor of the whitefly Bemisia tabaci interacts with the capsid proteins of begomoviruses and inhibits virus retention. Insect Mol Biol. (2023) 32:240–50. doi: 10.1111/imb.12827

PubMed Abstract | Crossref Full Text | Google Scholar

41. Aufauvre, J, Misme-Aucouturier, B, Viguès, B, Texier, C, Delbac, F, and Blot, N. Transcriptome analyses of the honeybee response to Nosema ceranae and insecticides. PLoS One. (2014) 9:e91686. doi: 10.1371/journal.pone.0091686

PubMed Abstract | Crossref Full Text | Google Scholar

42. Junell, A, Uvell, H, Davis, MM, Edlundh-Rose, E, Antonsson, A, Pick, L, et al. The POU transcription factor drifter/ventral veinless regulates expression of drosophila immune defense genes. Mol Cell Biol. (2010) 30:3672–84. doi: 10.1128/MCB.00223-10

PubMed Abstract | Crossref Full Text | Google Scholar

43. Liu, X, and Chu, Z. Genome-wide evolutionary characterization and analysis of bZIP transcription factors and their expression profiles in response to multiple abiotic stresses in Brachypodium distachyon. BMC Genomics. (2015) 16:1–15. doi: 10.1186/s12864-015-1457-9

Crossref Full Text | Google Scholar