Arachis hypogaea seedlings (cv. Georgia Green) were germinated and maintained in a growth room at 25–30°C, 40–50% relative humidity (rh) and L14:D10 photoperiod. Subsequently, sprouted seeds were transplanted and maintained in thrips-proof cages in a greenhouse at 25–30°C, 80–90% rh, and L14:D10 photoperiod. TSWV-infected volunteer A. hypogaea plants were collected from Belflower Farm (Coastal Plain Experimental Station, Tifton, GA, USA) and maintained in a greenhouse under the same environmental conditions. Leaflets with TSWV symptoms were collected, placed in a Munger cage, into which non-viruliferous F. fusca female adults were released (Munger, 1942). After the next generation adults emerged, potentially viruliferous adults obtained from the cages were used to inoculate prepared A. hypogaea plants. Ten female adults (up to 2 days old) were released on 1-week-old plants that had been dusted with 0.05 g of pine (Pinus taeda L.) pollen. Each plant was enclosed in a Mylar film (Grafix, Cleveland, PA, USA) cage with a copper mesh top. Plants were maintained in thrip-proof cages in the greenhouse as described above. Subsequently 3 weeks post thrips inoculation, infection status of the plants was determined by real-time quantitative reverse transcription polymerase chain reaction (qRT-PCR) using TSWV N-gene-specific forward and reverse primers, 5′-GCTTCCCACCC TTTGATTC-3′ and 5′-ATAGCCAAGACAACACTGATC-3′, respectively (Rotenberg et al., 2009). Plasmids with TSWV N-gene inserts were used as external copy number standards.

Adult F. fusca were initially collected from A. hypogaea blooms from Belflower farm. Subsequently, non-viruliferous and potentially viruliferous F. fusca colonies were separately reared on non-infected and TSWV-infected A. hypogaea leaflets, respectively, in Munger cages under the conditions described earlier. F. fusca reared for an entire generation (egg to adult) on TSWV-infected A. hypogaea leaflets alone were considered as potentially viruliferous. The TSWV infection status of subsamples from potentially viruliferous and non-viruliferous thrips used for this study was determined by qRT-PCR using the primers and conditions as described above. In both non-viruliferous and potentially viruliferous colonies, a majority of thrips (>95%) were females. Rotenberg et al. (2009) revealed that TSWV nucleocapsid RNA copies were two to three times greater in female thrips [Frankliniella occidentalis (Pergande)] when compared with male thrips. Although a different thrips species was the focus of this study [F. fusca (Hinds)], only females were used.

Approximately 160 female adults (up to 2 days old) were separately pooled from viruliferous and non-viruliferous thrips colonies, and total RNA was extracted. Similarly, leaf tissue from three TSWV-infected or non-infected plants (about 8 weeks old) was used for plant total RNA extraction. A total of two independent biological replicates were included for each F. fusca treatment and three independent biological replicates for each A. hypogaea treatment. Total RNA from thrips and plant samples were extracted using Direct-zol RNA Miniprep kit (Zymo Research, Irvine, CA, USA) using the manufacturer’s protocol. Library construction and sRNA sequencing (50 bp, single end) were carried out by the Beijing Genomic Institute using the Illumina HiSeq sequencing system. Samples were multiplexed and sequenced over two lanes. All sequence files are publicly accessible via the NCBI Sequence Read Archive (uninfected A. hypogaea: SRR3727343, SRR3727344, SRR3727345; TSWV-infected A. hypogaea: SRR3727346, SRR3727347, SRR3727348; uninfected F. fusca: SRR3727364, SRR3727365; TSWV-infected F. fusca: SRR3727366, SRR3727368).

Alignments of sRNA reads were carried out against the following TSWV RNA segment accessions: L RNA - NC_002052.1; M RNA - AY744483.1; S RNA - AY744475.1. For A. hypogaea sample alignments, reference sequences were obtained for A. ipaensis, a sequenced progenitor of A. hypogaea. These references included an annotated transcriptome (Araip.K30076.a1.M1) and mobile element data set (BB051914; Dash et al., 2016). Previously published raw sRNA sequence files from TSWV-infected S. lycopersicum (SRR947065) and N. benthamiana (SRR947064) were re-processed and used in this study for comparative purposes (Mitter et al., 2013). In this previous study, S. lycopersicum and N. benthamiana plants were mechanically inoculated with TSWV 3–4 weeks post-germination and leaf samples from TSWV positive plants were collected 17 days post inoculation.

Briefly, adapter sequences were removed from 50 bp single-end reads, which were then collapsed to identical sequences with maintained counts using the FASTX-Toolkit package. Only reads with lengths between 18 and 32 nt were retained for alignment normalization. Alignment of collapsed reads to reference sequences was performed using the SCRAM software package which allowed for exact matches to the reference sequence only (Fletcher, 2016). Normalization of aligned read count at each position was calculated as reads aligned per million reads between length 18 and 32 nt in the collapsed read file, allowing for robust quantitative comparison between alignments. Boxplots and sRNA length by count profiles were generated in R Studio using custom R scripts. Graphical representations of alignment profiles were produced either directly from the Python alignment software or via the Circos software package using tabulated output (Krzywinski et al., 2009). Scatter plots and Venn diagrams were generated in R Studio, and Euler diagrams in eulerAPE.

The size distributions of sRNAs in TSWV-infected and uninfected samples are shown in Figure 2A. The uninfected samples represent only the host sRNAs whereas the infected samples include total sRNAs (host and vsiRNAs). We have also aligned only vsiRNAs in the infected samples. The categorisation of 20–24 nt total sRNAs into host sRNAs (non-vsiRNAs) and vsiRNAs is presented in Supplementary Table S1. The bimodal size distributions of sRNA reads in healthy A. hypogaea and F. fusca are typical of their origin, with 21 and 24 nt sRNAs dominating in the plant host and 22 nt and 27–28 nt sRNAs most prevalent in the insect vector (Figure 2A). The relative abundance of 21 and 24 nt peaks apparent in the healthy A. hypogaea profile closely resembled that of N. benthamiana, a Solanaceae member, in a previous report (Mitter et al., 2013). Similarly, the distribution of sRNA reads in the infected and healthy F. fusca samples matched those in profiles of insect species as diverse as the silkworm Bombyx mori (L.; Jagadeeswaran et al., 2010) and migratory locust Locusta migratoria (Wei et al., 2009). L. striatellus infected with Rice stripe virus also displayed a bipartite profile, though a peak was evident at 26–27 nt rather than 27–28 nt, and contrastingly, the abundance of the sRNAs making up the larger sized peak was lower than that of the 22 nt sRNAs (Xu et al., 2012). The distinct profiles of the plant host and insect vector likely reflect the differences in sRNA biogenesis and function among the organisms. The use of female thrips may have influenced the prevalence of vsiRNAs observed, as they could be more abundant in female thrips than in male thrips (Rotenberg et al., 2009).

Remarkably, the impact of TSWV infection on the size distribution of total sRNAs in A. hypogaea is profound, with a significant reduction in total 24 nt sRNAs after infection (Figures 2A,B). Compared to uninfected plants, there is also an increase in 22 nt sRNAs, and to a lesser extent, 20 nt sRNAs (Figure 2B). Such a change in relative abundance of 22and 24 nt sRNAs due to TSWV infection is in contrast to an increase in 24 nt sRNAs and decrease in 22 nt sRNAs in other reported species, such as Bamboo mosaic virus-inoculated N. benthamiana plants (Lin et al., 2010). A slight decrease in 24 nt sRNAs in cotton plants infected with Cotton leaf roll dwarf virus indicates a similar pattern, but at a much reduced magnitude compared to TSWV-infected A. hypogaea (Silva et al., 2011). These data suggest that infection by TSWV significantly alters the productive balance of the plant’s sRNA biogenesis machinery. In contrast, infection of the F. fusca vector by TSWV results in minimal change in the relative abundance of the various sRNA discrete sizes (Figure 2A).

Due to their role in maintenance of genome stability via repression of transposable elements, the reduction in the number of 24 nt sRNAs in TSWV-infected A. hypogaea could have significant physiological consequences. Additionally, an increase in 22 nt sRNAs may lead to other further changes in endogenous sRNA abundance, as these sRNAs trigger biogenesis of phased sRNAs spanning complementary RNA sequences, along with silencing of the affected genes (McHale et al., 2013). Widespread perturbation of gene expression in infected plants, either directly or indirectly, may be a fundamental factor in the infected plant phenotype.

Analysis of sRNA alignments to repetitive and transposable elements within the A. hypogaea genome demonstrated a clear reduction in sRNAs derived from these elements in TSWV-infected plants (Figure 3). This shift likely has wide-ranging impacts, with transposition events potentially happening at a greater frequency, and the expression of nearby endogenous genes altered. In support of this prediction, previous studies have shown that the disruption of biogenesis of 24 nt sRNAs in plants had an extensive impact on the degree of DNA methylation along with the activity of targeted transposable elements, and the expression of neighboring genes (Jia et al., 2009; Groszmann et al., 2011).

Alignments to protein-coding sequences also demonstrated clear changes in the abundance of 21, 22, and 24 nt sRNAs targeting these sequences (Figure 4). Interestingly, an increase in the variance of coverage of coding sequences by 21, 22, and 24 nt sRNAs is particularly evident, pointing toward possible broad up- and down-regulation of gene expression in infected leaf tissue. A previous study has suggested that the TSWV-derived Non-structural protein (NSs) binds directly to intermediates of sRNA processing and potentially suppresses loading of sRNA into AGO complexes required for degradation or translational inhibition of complementary RNAs (Schnettler et al., 2010). Therefore, the disruption of sRNA levels observed in this study and the potential for NSs action on plant-generated sRNAs may be a means for perturbation of target RNA transcript levels and consequential physiological and developmental changes.

This report adds further evidence for the species-specific impact virus inoculation has on the dynamics of sRNA processing. The exact mechanism driving this shift in sRNA abundance and the reason it differs between species remains to be fully elucidated. DCL2, DCL3, and DCL4 homologs of plants all cleave dsRNA substrates, either endogenous or exogenous/viral derived dsRNAs, into 22, 24, and 21 nt sRNAs, respectively (Henderson et al., 2006; Bologna and Voinnet, 2014). Studies on Arabidopsis dcl mutants show DCL4 and DCL2 are partly redundant and act hierarchically in processing dsRNA substrates associated with post-transcriptional gene silencing (Parent et al., 2015). In contrast, DCL3 is accepted as the primary enzyme for generation of 24 nt siRNAs from dsRNAs produced by the Pol IV-RDR2 complex and associated with transcriptional silencing of repetitive elements and transposons (Gasciolli et al., 2005; Pontes et al., 2006; Wierzbicki et al., 2009; Johnson et al., 2014). Our data demonstrated that TSWV infection in A. hypogaea, but not the two Solanaceae species (Mitter et al., 2013), triggered a genome-wide depletion in 24 nt siRNAs associated with repetitive elements and transposons (Figure 3). This suggests that in A. hypogaea there was a TSWV-dependent suppression in expression levels, subcellular localization and/or the enzyme function of DCL3, or another upstream component of the transcriptional gene silencing pathway, such as Pol IV or RDR2 (Wierzbicki et al., 2009). Thus, inhibition of DCL3 or another component of the transcriptional gene silencing pathway is the most likely explanation for the shift in production of 24 nt siRNAs to 22 nt siRNAs in TSWV-infected A. hypogaea. Such a reaction to TSWV infection on the accumulation of 24 nt siRNAs clearly varies between plant species. Further expression analysis on components of the transcriptional gene silencing pathway may shed light on differential effects on TSWV infection in different host species.

The portion of total sRNAs that align to the infecting viral genome has been shown to vary greatly depending on the virus and the plant host (Dehaan et al., 1991; Bolduc et al., 2010; Mitter et al., 2013; Hong et al., 2015). Approximately 16% of all reads in each infected A. hypogaea replicate aligned to the TSWV genome, while just 0.4% aligned to the virus in the infected F. fusca vector (Table 1). The abundance of vsiRNAs aligning to the TSWV genome in A. hypogaea is high compared to S. lycopersicum and N. benthamiana infected with the same viral isolate (6% and 0.02%, respectively), though in the range of other host–virus combinations such as peach infected with Peach latent mosaic viroid (11.2 %; Bolduc et al., 2010) and rice infected with Rice dwarf phytoreovirus (∼22%; Hong et al., 2015).

Tospoviruses are characterized by a tripartite genome structure consisting of large (L), medium (M) and small (S) single-stranded RNA segments with sense and ambi-sense coding frames. The L segment consists of a large open reading frame (ORF) that encodes RdRp (RNA-dependent RNA polymerase), which is required for viral replication (Dehaan et al., 1991). ORFs for Gn/Gc (precursor for the N-terminal G2 and C-terminal G1proteins) and NSm (non-structural protein) occur in the M segment. NSm is required for cell-to-cell viral movement (Soellick et al., 2000), and Gn/Gc for synthesis of glycoproteins present in the viruses’ lipid envelope (Snippe et al., 2007). The S segment encodes the N (nucleocapsid) gene, which encapsidates the RNA segments (Dehaan et al., 1990), and NSs (non-structural/suppressor of silencing protein), which is of particular importance in the defense against the plant’s RNAi defense strategy. Indeed, the degree of NSs protein accumulation in plant tissues has been associated with the severity of TSWV infections (Takeda et al., 2002). NSs can bind to both dsRNA and vsiRNA targets, allowing it to block antiviral RNA silencing before and after DCL-mediated dsRNA cleavage (Schnettler et al., 2010). Importantly, NSs has been shown to be effective in not only plant hosts, but also in insect vectors (Margaria et al., 2014), possibly through its affinity for longer dsRNA forms. Alignment of reads from infected plant samples to individual genome segments of TSWV has been previously reported, also with dissimilar results depending on the plant species and viral isolate (Mitter et al., 2013; Margaria et al., 2015).

In this study, TSWV-aligned vsiRNAs from A. hypogaea demonstrated a distribution distinct from vsiRNAs from the Solanaceae hosts, with a generally high portion aligned to the L segment and variable percentages aligning to the M and S segments, depending on the replicate. In contrast, F. fusca-derived vsiRNAs more closely matched the aligned vsiRNA distribution from S. lycopersicum and N. benthamiana (Figure 5; Table 2).

The use of multiple fully independent biological replicates for A. hypogaea sRNA sequencing provides a detailed picture of the strand polarity of alignments to individual TSWV gene sequences (Figure 6). Significant variation amongst the A. hypogaea replicates indicates that there is little preference for a particular alignment orientation. This is in contrast to earlier reports for TSWV, where a greater abundance of vsiRNAs generally aligned to the sense strand of the virus in Solanaceae hosts (Mitter et al., 2013). The possibility exists that greater replication of sequencing for other TSWV-infected plant species may demonstrate similar variation. F. fusca, on the other hand, shows a sense preference for the Gn/Gc, RdRp and Nc genes, and an antisense preference for the NSm and NSs genes.

Due to vsiRNA being processed from dsRNA intermediates, it could be expected that vsiRNAs are equally distributed toward sense and antisense strands. The observed strand bias of vsiRNAs in F. fusca may be due to vsiRNAs binding to highly abundant complimentary viral transcripts, as has been suggested in previous studies (Smith et al., 2010; Harris et al., 2015). Such binding could potentially sequester and remove vsiRNA from the pool for sRNAs submitted for sequencing. If such an explanation were correct, it is predicted that the difference between strand biases observed for A. hypogea compared to F. fusca may be due to host-specific differences in the level of viral RNA transcript expression and/or availability.

The location and magnitude of vsiRNA hotspots spanning the TSWV genome in A. hypogaea and F. fusca samples are partially similar, though key differences exist (Figure 7A). In general, 21 and 22 nt vsiRNA hotspots within each species are co-located, with 21 nt peaks of greater magnitude in A. hypogaea, and 22 nt peaks larger in F. fusca. Visualization of vsiRNA coverage via circos plots demonstrates that relative to F. fusca, hotspots are more evenly dispersed across each genome segment in A. hypogaea. The F. fusca plots show fewer easily discernible peaks in the large RNA segment, and a greater abundance in the sRNA segment. Hotspots often occur in coding regions of the genome for both species, with a distinct gap apparent between adjacent genes located on the sRNA segment. Whilst hotspots do sometimes occur in the same locus in both the plant host and insect vector, many are unique to each organism; reinforcing that RNAi targeting mechanisms do have species-specific preferences. Importantly, these differences also indicate that TSWV-specific sRNAs in the insect are not simply derived from the plant itself via feeding, and are unique to the insect vector.

A comparison of A. hypogaea vsiRNA hotspots with those previously reported for S. lycopersicum and N. benthamiana infected with the same TSWV isolate (Mitter et al., 2013) shows distinct differences between plant families (Figure 7B). Similar to F. fusca, the vsiRNA hotspots in the L RNA segments corresponding to the RdRP gene are lacking in S. lycopersicum and N. benthamiana, indicating a very specific focus for A. hypogaea. Accordingly, though the external symptoms of TSWV infection are comparable in severity between A. hypogaea and S. lycopersicum, this A. hypogaea-specific direction of vsiRNAs toward the TSWV’s RdRP sequence may have significant downstream impacts upon pathogenicity relative to the Solanaceae.

The biogenesis and specific effects of vsiRNAs that form hotspots in alignment profiles are an area of ongoing study. Highly expressed vsiRNAs are thought to be a result of advantageous structural conditions in the confined regions of dsRNAs that vsiRNAs are generated from (Zhang et al., 2015). Differences in vsiRNA hotspot prevalence among plant families infected with the same viral isolate indicate that the localized environment the virus inhabits varies enough to impact RNA secondary structure, and thus the relative makeup of vsiRNAs generated.

The presence and absence of discrete vsiRNAs that are unique or shared between host and vector, or between different plant species infected with TSWV, reinforce the interpretations of the hotspot analyses. Unique vsiRNAs are readily apparent in F. fusca samples, further supporting the indication that vsiRNA generation takes place within the thrips, and not solely through uptake via feeding on vsiRNAs produced in infected A. hypogaea leaf tissue (Figure 8A). Additionally, the differences between A. hypogaea and the Solanaceae remain pronounced (Figure 8B). Due to the low TSWV vsiRNA load in N. benthamiana, almost all vsiRNAs present are common to S. lycopersicum. However, a similar number of discrete vsiRNAs only show partial overlap between A. hypogaea and N. benthamiana/S. lycopersicum, reinforcing the idea of a distinctive RNAi-mediated response to TSWV infection for different host plant families.