Characteristics of vsiRNAs
To determine the roles of vsiRNAs in the interactions between the virus and the host plant, high-throughput sequencing of small RNAs was performed on BGISEQ-500 to analyze the vsiRNA population in BYDV-GAV-infected Xiaoyan 6 wheat plants. A total of 26,286,625 clean reads were yielded from the library for further analysis. Next, we mapped the clean reads to the viral dsRNA genome with 2 mismatches and obtained 2,440,544 filtered reads, which accounted for 9.28%. After processing and redundancy screening, 113,884 filtered siRNAs were obtained, which were recognized as vsiRNAs. The length distribution of vsiRNAs indicated that 21- and 22-nt vsiRNAs were abundantly present followed by 24-nt vsiRNAs in the library (Fig. 2a). The results showed that BYDV-GAV infection produced a significantly increased population of 21- and 22-nt siRNAs, suggesting that DCL4 and DCL2 enzymes play major roles in the biogenesis of BYDV-GAV-derived siRNAs.
Previous studies reported that the specific recruitment of sRNAs by distinct AGO proteins largely depends on the 5’-terminal nucleotide [50-51]. AGO2 and AGO4 primarily recruit small RNAs with a 5'-terminal A, whereas 5'-termini with U and C are preferentially recruited by AGO1 [52]. To investigate which AGO is primarily involved in the processing of BYDV-GAV-derived vsiRNAs, the 5’-terminal nucleotides of the vsiRNAs were analyzed. For the heavily accumulated 21- and 24-nt lengths, A and C accounted for a dominant proportion, while U and G were less abundant. According to the results, the BYDV-GAV-derived vsiRNAs had a 5’-terminal nucleotide bias toward A and U, which might be preferentially recruited by wheat AGO1, AGO2 and AGO4 to form the RISC (Fig. 2b).
To elucidate the origin of the BYDV-GAV-derived vsiRNAs, the polarity distribution of vsiRNAs in the positive (+) and negative (-) strands of the viral dsRNA genome was characterized by aligning the vsiRNAs with the genome. The results indicated that the hotspots of sense- and antisense-vsiRNAs that were derived from positive-sense (+)- and negative-sense (−)-strands of BYDV-GAV genomic RNA were nearly the same (Fig. 2c). However, the number of sense-vsiRNAs was slightly greater than that of antisense-vsiRNAs.
From the distribution of vsiRNA in positive-sense (+)- and negative-sense (−)-strands of BYDV-GAV genomic RNA, we observed that vsiRNAs were not continuously distributed through the genomic segments. Among all segments, vsiRNAs predominantly accumulated in ORF1, OPF3, and ORF4, which were primarily responsible for viral replication, infection, and movement, suggesting that vsiRNAs play a vital role in the processing of viral invasion.
Prediction of putative target transcripts of vsiRNAs
Target identification is beneficial to provide an in-depth understanding of the functions of vsiRNAs. In this study, vsiRNAs were employed to predict putative wheat target genes by the target analysis server TargetFinder with default parameters. A total of 23,719 target transcripts of 11,384 vsiRNAs were acquired after deleting duplicate items, of which 4,324 vsiRNAs had corresponding targets in wheat plants (Additional file 2: Table S5). Annotation of the predicted target genes was performed referring to the wheat transcript database from ensemble plants.
GO and KEGG pathway enrichment analysis of predicted targets
To determine the roles of the predicted target genes in wheat plants, GO and KEGG pathway enrichment analyses were performed to investigate the underlying functions. The GO enrichment analysis demonstrated that cellular and metabolic processes were the most significantly represented groups in the biological process category, cells, and cell parts were the most commonly enriched terms in the cellular component category, while catalytic activity and binding were the predominantly enriched groups within the molecular function category (Fig. 3a). The KEGG pathway enrichment analysis showed that the target genes were most abundant in enriched plant-pathogen interaction pathways (Fig. 3b). From the KEGG pathway, we also observed that the porphyrin and chlorophyll metabolism pathway was enriched, containing approximately 178 genes. We speculate that this pathway was involved in the yellowing of BYDV-GAV-infected wheat.
Expression of genes related to yellowing of BYDV-GAV-infected wheat
To identify the genes related to symptoms of yellowing, we selected a series of predicted targets involved in chlorophyll biosynthesis. Next, RT-qPCR was conducted to measure their relative expression level in wheat samples at 14 d after infection with BYDV-GAV, and we found that most of the targets showed a significantly decreased expression level. For example, the relative accumulation levels of cytochrome P450, chlorophyll synthase, chlorophyll a-b binding protein, glutamyl-tRNA reductase, glucose-6-phosphate isomerase, and chloroplast Mg-chelatase subunit were decreased by approximately 0.17-, 0.1-, 0.08-, 0.21-, 0.19-, and 0.14-fold, respectively (Fig. 4). The results indicated that infection with BYDV-GAV can inhibit the expression of some chlorophyll biosynthesis-related genes in wheat plants. Therefore, further study is urgently needed to confirm whether vsiRNAs had a downregulation effect on these predicted target genes. All of these predicted target mRNA cleavage sites of the vsiRNAs are listed in Additional file 3 Table S6.
ChlG-silenced wheat displays leaf yellowing
Among the targets, ChlG was predicted to encode chlorophyll synthase, which functions in the last step of the chlorophyll biosynthetic pathway, playing a decisive role in chlorophyll content. To determine the functions of the target gene, transient silencing of ChlG was implemented through virus-induced gene silencing (VIGS) with the BSMV system. At 10 d postinfiltration (dpi), we determined that the plants infiltrated by γ-PDS exhibited significant photobleaching, whereas the newly emerged leaves of ChlG-silenced wheat plants were yellowing, which is the typical symptom of BYDV-GAV (Fig. 5a). Meanwhile, the expression levels of ChlG and PDS genes in silenced plants were detected by RT-qPCR. The results showed that the expression level of the PDS gene was decreased by almost 75% and that the ChlG gene was reduced by 87% when compared to the negative control, which was inoculated with the empty vector (Fig. 5b).
After 7 d of infection with viruliferous aphids carrying BYDV-GAV, total RNA of the inoculated leaves was extracted to analyze the relative content of the virus. The CP gene of BYDV-GAV was employed to determine viral accumulation. The RT-qPCR results demonstrated that the accumulation levels of BYDV-GAV in ChlG-silenced wheat plants showed a greater than 4.5-fold increase compared to the negative control plants (Fig. 5c). The increased viral accumulation in ChlG-silenced wheat plants indicated that ChlG might be involved in the resistance to BYDV-GAV infection.
Alteration of chlorophyll content in wheat and ChlG-silenced wheat plants
To further confirm the relationship between chlorophyll content and the expression of the ChlG gene in wheat plants, we measured the chlorophyll content by utilizing the absorbance at 663.6 and 646.6 nm to calculate the final quantitative value. Both wheat and ChlG-silenced wheat plants infected with viruliferous aphids were collected to analyze the chlorophyll content. The chlorophyll content in BYDV-GAV inoculated wheat plants showed a significant decrease of approximately 32% compared with mock-inoculated plants at 10 dpi. As expected, this accumulation was also significantly reduced by 38% in the ChlG-silenced wheat plants compared to the mock-inoculated plants at 10 dpi (Fig. 6). The results indicated that the yellowing symptom of BYDV-GAV might be closely associated with the destruction of ChlG.
vsiRNA8856 from the common region of MP and CP
Previous studies have reported that vsiRNAs can target host endogenous gene transcripts in plants [20, 21, 27]. To explore the relationship between the decreased expression of ChlG and the BYDV-GAV-derived vsiRNAs, we analyzed the results of the vsiRNA that was predicted to target ChlG. We found the cleavage site primarily at the 30-, 75-, and 1097-bp sites of the ChlG transcript. Then, the expression of vsiRNAs was detected in wheat plants by RT-PCR (Fig. 7). In addition, we found the production of large amounts of vsiRNA8856 along with BYDV-GAV infection, and 346 reads were detected from the sequencing data.
Validation of vsiRNA target by 5’ RLM-RACE
Although our results show that downregulation of ChlG correlated strongly with symptoms of BYDV-GAV, it is not known how viral infection downregulates host genes. To identify the accuracy of predicted vsiRNA targets, 5’ RLM-RACE was used to confirm the interaction between vsiRNAs and predicted target genes. Then, cleavage sites were found between the 10th and 11th nucleotides of vsiRNA8856 for the ChlG transcript. The results of 5’ RLM-RACE further indicate that the vsiRNAs cleave their corresponding target genes at specific sites (Figure 8), which suggested that the downregulation of ChlG contributed to cleavage by vsiRNAs from BYDV-GAV.
Cloning, multiple alignments and phylogenetic analysis of ChlGs
The full-length ChlG gene was cloned from wheat plants by RT-PCR. The results showed that the sequence of the ChlG gene is 1,134 bp with eight ORFs, of which ORF1 contains the complete sequence. The ChlG protein contains a total of 377 amino acids (aa) with a molecular weight of 40.6 kDa and a theoretical pI of 7.59. The subcellular localization prediction suggests that the ChlG gene may be located at the plasma membrane.
To further evaluate the phylogenetic relationships between wheat ChlG and other homologous species, we aligned the aa sequences among the species Hordeum vulgare (HvChlG), Oryza sativa (OsChlG), Brachypodium distachyon (BdChlG), Sorghum bicolor (SbChlG), Setaria italica (SiChlG), Zea mays (ZmChlG), Avena sativa (AsChlG), Glycine max (GmChlG) and Triticum aestivum (TaChlG) (Figure 9a). The sequence alignment results showed that ChlG shares high aa sequence identities. In addition, we constructed a phylogenetic tree using these aa sequences (Figure 9b). The results indicated that TaChlG is most closely related to HvChlG followed by BdChlG and AsChlG. It has been indicted that BYDV-GAV can infect Hordeum vulgare, Avena sativa, and Brachypodium distachyon and cause symptoms. Hence, we infer that ChlGs have similar functions that were targeted by BYDV-GAV for the production of symptoms in these closely related species.