Identification of Novel and Conserved microRNAs in Homalodisca vitripennis, the Glassy-Winged Sharpshooter by Expression Profiling

The glassy-winged sharpshooter (GWSS) Homalodisca vitripennis (Hemiptera: Cicadellidae), is a xylem-feeding leafhopper and an important vector of the bacterium Xylella fastidiosa; the causal agent of Pierce’s disease of grapevines. MicroRNAs are a class of small RNAs that play an important role in the functional development of various organisms including insects. In H. vitripennis, we identified microRNAs using high-throughput deep sequencing of adults followed by computational and manual annotation. A total of 14 novel microRNAs that are not found in the miRBase were identified from adult H. vitripennis. Conserved microRNAs were also found in our datasets. By comparison to our previously determined transcriptome sequence of H. vitripennis, we identified the potential targets of the microRNAs in the transcriptome. This microRNA profile information not only provides a more nuanced understanding of the biological and physiological mechanisms that govern gene expression in H. vitripennis, but may also lead to the identification of novel mechanisms for biorationally designed management strategies through the use of microRNAs.


Introduction
MicroRNAs (miRNAs) are a class of non-protein coding small RNAs that are known to play important roles in post-transcriptional level gene regulation and expression (reviewed in [1,2]. MicroRNAs are typically 19-25 nt in length and are the most abundant class of endogenous small RNAs in both plants and animals [3]. MicroRNAs can also originate from either noncoding or coding regions of transcripts. The primary transcripts of microRNAs are predominantly cleaved and processed within the cell nucleus by the nuclease Drosha to generate precursor microRNAs (pre-miRNAs) with a characteristic hairpin like secondary structure. The pre-miRNAs are then transported into the cytoplasm by Exportin-5, where the terminal loop regulation; and to identify possible targets for RNAi-based insect control. MicroRNA analysis of H. vitripennis can also provide insights into potential roles of microRNAs in insect biology and pathogen-host interactions; that can generally aid us in understanding the molecular mechanisms of miRNA regulation.
The objectives of this study were to use a combination of deep sequencing of small RNAs of H. vitripennis and computational analysis to systematically predict microRNA species in the H. vitripennis and to identify their possible targets. We performed Illumina-based sequencing of small RNAs isolated from adult H. vitripennis and developed better understanding of the complexity of microRNA regulation in H. vitripennis. Our analyses identified a number of new and novel microRNAs and their precursor RNAs. Our results provide a global view of the micro-RNAome of adult H. vitripennis and pave the way for its further analysis and exploitation.

Materials and Methods
Maintenance of H. vitripennis and total RNA extraction A colony of H. vitripennis was maintained at the University of California-Davis Contained Research Facility (CRF) in cages containing a mixture of host plants as previously described [31]. Total RNAs were extracted using TRIzol reagent (Invitrogen, Carlsbad, CA, USA) from eight day-old adults (males and females) or 2 nd -3 rd instar individual nymphs Library preparation and sequence data generation Total RNA (2 μg per library) was used as the template to construct single-end indexed Illumina small RNA libraries. Samples were processed according to Illumina's TruSeq Small RNA sample preparation guide. Briefly, the small RNAs were ligated to Illumina's small RNA 3' and 5' adaptors. Subsequently cDNA was synthesized by reverse transcription and amplified by PCR (11 cycles) purified by two size selection gels. Sequencing of the small RNA library was performed by the Illumina Genome Analyzer. The small RNA raw sequence data in the form of a FASTQ file has been uploaded to the SRA database at NCBI (accession number: SAMN03853644; Accession ID: SRS985993) (http://www.ncbi.nlm.nih.gov/biosample/?term= SAMN03853644). The reference numbers for the SRA experiments are SRX1090328 and SRR2095935 (http://www.ncbi.nlm.nih.gov/sra/SRR2095935/). The Bio Project accession number allotted to the sequence information is PRJNA289427 (Accession: SRP060579).

Data extraction and analysis
All of the Illumina sequencing data were initially converted to FASTA format from FASTQ. Small RNA sequences were extracted from raw reads by removal of the adaptor sequences. Based on the length of the mature miRNA and adapter length, small RNA sequences of 22 nt in length were filtered out, checked for redundancy and were queried against the ribosomal and transfer RNA databases. The sequences were then queried against the H. vitripennis transcriptome [15] and the sequences that mapped to the transcriptome were retained. The sequences were also aligned against miRBase [16,17,[33][34][35][36] to identify the conserved (two mismatches allowed) and novel microRNAs. Mature miRNA sequences were used as query against H. vitripennis transcriptome to locate corresponding precursors that had at least 18 matched base pairs, one central loop and folding energy lower than 18 kCal/mol. The frequency of the miRNAs and miRNA Ã s were calculated based on the read numbers in the library. All bioinformatics analyses were performed using custom written perl scripts and through the help of the scripts made available through the miRDeep2 software [11]. Further annotation of the micro-RNA sequences was performed through manual annotation. We were able to identify precursors for the newly identified microRNAs from the transcriptome information of H. vitripennis [15]. This is important because we did not have to use a proxy reference as suggested in the recent literature on the correct usage of microRNA discovery pipelines in non-model organisms [37].

MicroRNA target predictions
Target prediction for animal miRNAs is complex because of imperfect complementarity between miRNAs and their mRNA targets. The target candidates for miRNAs of H. vitripennis were chosen based on our mRNA dataset (~52,700 sequences) that were generated through transcriptome sequencing [15]. The miRNA target candidates were analyzed with target prediction software MiRanda version 3.3a [38]. The three criteria that were allowed for the screening of miRNA target alignments were: no mismatch at the seed region (positions 2 to 7 from the 5' end of predicted miRNA); not more than 1 G:U paring at the seed region; and not more than 1 gap in the miRNA:mRNA duplex. Except for the energy and score parameters, the remaining parameters in the MiRanda version 3.3a software were set to default settings. The total score cutoff value was set for !145 and energy threshold value of −10 kcal/mol was set as second parameter.

Quantitative real-time PCR
Total RNA (isolated as described before) from adult and nymphal (2 nd -3 rd instars) H. vitripennis was used to study the expression pattern of the conserved and novel microRNAs. The miRNA expression was measured using a two-step process. In the first step, a stem-loop (RT) primer (designed based on previous reports [39]) was hybridized to the miRNA and reverse transcribed in a pulsed RT reaction [40]. In the second step, the RT reaction product was PCR amplified using a miRNA specific forward primer and a universal reverse primer (S1 Table) in real time with SYBR green chemistry using a Bio-Rad CFX Real-Time PCR detection system [40]. Quantification of the relative changes in miRNA expression was performed using the method previously described [40]. In this relative quantification method, the sample reference was chosen as miR10a and the endogenous control was chosen as ubiquitin [41]. The data for relative quantities were converted to fold differences by logarithmic transformation to express the data as a normal distribution. The data presented are the averages of three measurements.

MicroRNA northern blots
Some of the conserved and novel microRNAs were validated using microRNA northern blots through the use of 5' end labeled probes. Total RNA was isolated from the individual adult H. vitripennis using Trizol reagent and small RNAs were precipitated using 5 M NaCl and 50% PEG8000 [42]. Briefly, three micro grams of adult H. vitripennis small RNAs were separated on denaturing 8 M urea-15% polyacrylamide gels and transferred to a nylon membrane (Amersham Hybond™-NX) as described previously [43]. The membrane was then UV cross-linked (Stratagene). Synthetic antisense DNA oligomers were 5'-end labeled using [γ 32 P]-ATP (Perkin Elmer) and T4 Polynucleotide kinase (New England Biolabs) and used as probes for hybridization. The probes were purified on a Sephadex G-25 column (mini quick spin oligo columns, Roche Diagnostics). Prehybridization was carried out for 5 hours using ULTRAhyb-oligo hybridization buffer (Ambion, AM8663). The hybridization was performed using a probe specific activity of 2 x 10 6 cpm per ml of hybridization buffer. The hybridization was carried out for 16 hours at 42°C in a hybridization oven. Following hybridization, the membrane was washed twice by shaking in low stringency wash solution (2X SSC, 0.1% (w/v) SDS) at 42°C for 20 min. The blots were removed, drained, wrapped in Saran Wrap and exposed to X-ray film at −80°C for 6-8 days. MicroRNA Marker (New England Biolabs, Ipswich, MA) was used on the same gel to estimate RNA sizes. Wild type potato (Solanum tuberosum cv. Desiree) small RNA samples on which the insects fed were used as negative controls.

Overview of the dataset
Sequencing of the adult H. vitripennis small RNA library generated 22 million reads. After small RNA read processing, the remaining sequence reads (43% of the total reads) were mapped to an artificial build of the H. vitripennis transcriptome [15]. Available scripts from the miRDeep2 software package were used to identify novel microRNAs [11]. The miRDeep2 analysis identified microRNA sequences along with their abundance. We were also able to investigate the secondary structure of each potential precursor microRNA as well as the coordinates of those precursor microRNAs in the transcriptome of H. vitripennis. In the filtering step, potential precursors that were grossly inconsistent with miRNA biogenesis were discarded. Based on the algorithmic scores assigned from the scripts, each potential microRNA precursor was verified for the combined compatibility of energetic stability, position and frequency of reads with Dicer processing.

Novel microRNAs in H. vitripennis
Novel species-specific microRNAs can be identified and annotated on the basis of 1) criteria defined for the identification of novel microRNAs (i.e., the presence of a miRNA-star strand, low-energy, and fold back precursor structure [30,44]), and 2) the lack of resemblance to known microRNAs in the available miRBase database. We used scripts from the miRDeep2 software package [11] to identify novel microRNAs in H. vitripennis. Because A. pisum is the closest hemipteran species (based on the available sequence information) to H. vitripennis, microRNAs from A. pisum were used as reference sequences. Additional microRNA sequences from B. mori, T. castaneum, and N. vitripennis were used for more accurate guidance. The precursor sequences and structures were identified from the available transcriptome information of H. vitripennis. Manual analysis of the algorithmic predictions helped us to identify and annotate 14 novel microRNAs ( Table 1). All of these sequences have low-energy, fold-back precursor structures. The consensus star strand and the consensus precursor  sequences for the microRNAs were identified (Table 2). Interestingly, a majority (11 out of 14) of the novel microRNAs that were identified in H. vitripennis showed a 5' terminal uridine bias. This 5' terminal uridine bias for microRNAs was previously reported to allow for specific Argonaute 1 (AGO1) recognition and loading in the RISC complex. Two of the novel micro-RNAs, Hvi-miR2062 and Hvi-miR6225, had 5' terminal guanosine and one Hvi-miR24657 had a 5' terminal adenosine. Recent studies on miRNA loading show an increase in loading of microRNAs with a 5' terminal adenosine bias to AGO2 and those with a cytidine bias into AGO5 [6].
Survey of the potential microRNA encoding genes in H. vitripennis A total of 372 mature miRNAs from A. pisum, B. mori, T. castaneum, and N. vitripennis were used as reference sequences to identify potential microRNA encoding genes in H. vitripennis.
The conserved and newly identified novel microRNA sequences were used to search the transcriptome assembly of H. vitripennis using scripts from bowtie [45] and mirdeep2 [11] functional algorithms, and through UNIX-based scripts. Under these criteria, a total of 14 transcripts that potentially encoded the newly identified microRNAs were identified from the H. vitripennis transcriptome ( Table 3). The sequence from each locus was used to determine the existence of a stem-loop pre-miRNA structure using mfold [46]. The microRNA target loci annotations from the transcriptome data were used to identify the relevance of the potential origins of the microRNAs in the transcriptome. Based on this analysis, seven (Hvi-miR13059, Hvi-miR9237, Hvi-miR19117, Hvi-miR24402, Hvi-miR41359, Hvi-miR91 and Hvi-miR24657) of the 14 identified microRNAs have origins that correspond to transposable elements (Table 3). Of these seven, Hvi-miR24402 and Hvi-miR9237 originated from loci that belonged to the TC3 type of transposable elements ( Table 3). One of the microRNAs that did not originate from a transposable element sequence, Hvi-miR2062, originated from a gamma amino butyric acid (GABA) receptor gene locus [47]. The coordinates of the originating loci of the precursor microRNAs in the transcriptome of H. vitripennis are presented in Table 3. Energy and fold-back structures of these precursor microRNA sequences, containing the mature strand and the consensus star strand sequences, are shown in Table 4.

Expression of microRNAs in H. vitripennis
Based on our analysis, we were able to identify conserved and novel microRNAs from adult H. vitripennis (Fig 1 and Table 1). Further validation of their expression was done using stemloop quantitative RT-PCR and microRNA northern blots. We tested six of the different micro-RNAs using the primer sets described in S1 Table and Fig 2. We were able to validate the expression of six of the conserved microRNAs from H. vitripennis (Hvi-miR1692, Hvi-miR159, Hvi-miR10a, Hvi-miR184, Hvi-miR263, Hvi-miR3256) and six novel microRNAs (Hvi-miR28196, Hvi-miR9237, Hvi-miR41359, Hvi-miR91, Hvi-miR13509 and Hvi-miR24402) through stem-loop quantitative RT-PCR (Fig 2). We were also able to detect the microRNAs in nymphal (2 nd and 3 rd instars) H. vitripennis (Fig 2). The expression of Hvi-miR171, Hvi-miR276, and Hvi-miR29828 was validated by northern hybridization of small RNAs isolated from individual adult H. vitripennis (Fig 3). The expression of the conserved microRNAs (Hvi-miR171, Hvi-miR276, and Hvi-miR29828) was confirmed in all of the adults tested (5 or 6 individuals). The northern blot results validated the expression of these conserved and novel microRNAs in adult H. vitripennis (Fig 3).

Transcriptome wide survey of the potential microRNA targets in H. vitripennis
MicroRNAs are known to regulate the mRNA levels thus influencing the development and immunity of the host organism [48][49][50][51]. The majority of microRNAs in animals are shown to interact with their targets through their 3' and 5'-untranslated region (UTRs). In order to further understand the biological functions of novel microRNAs in adult H. vitripennis, we screened the transcriptome of H. vitripennis [15] for the potential miRNA targets (~52,700 putative mRNAs) using the microRNA target prediction program Miranda 3.3 with threshold parameters as described in materials and methods [38]. On average, 9,132 potential targets of the 14 novel microRNAs were identified (Table 6). There were on average of 608 targets per novel microRNA ( Table 6). The number of potential targets varied in range from 172 (Hvi-miR66) to 1,583 (Hvi-miR13059) ( Table 6). Each microRNA was predicted to target multiple mRNAs; the top mRNA targets for each of the novel microRNAs along with their description and scores are shown in Table 7. The predicted targets of the novel microRNAs were involved in a diverse range of functions that include developmental processes, reproduction, enzymatic process, and ion channels (Table 7 and S2 Table).

Discussion
MicroRNAs form an important class of small RNAs that play important roles in insect development [8]. Here we identified 14 novel microRNAs from adult H. vitripennis. The precursors of these novel microRNAs were identified using previously determined and publicly available transcriptome information for H. vitripennis. These precursors were found to form low-energy, fold-back structures. The concise coordinates of the identified precursor microRNAs helped us to excise the precursor sequences for structural analysis. Furthermore, we also identified conserved microRNAs from H. vitripennis by homology analysis with microRNAs from other insects in miRBase (Fig 1). Sequence analysis of the 14 novel microRNAs from H. vitripennis suggests a strong bias for a 5' terminal uridine. This 5' terminal bias was found in 11 of the 14 novel microRNAs that we identified (Hvi-miR29035, Hvi-miR13059, Hvi-miR28196, Hvi-miR9237, Hvi-miR19117, Hvi-miR24402, Hvi-miR66, Hvi-miR41359, Hvi-miR29828, Hvi-miR20768 and Hvi-miR91). In two other microRNAs (Hvi-miR2062 and Hvi-miR6225) we observed a guanosine bias whereas Hvi-miR24657 showed a bias for adenosine. This difference in the 5' terminal nucleotide leads to a loading bias of the microRNA into different Argonaute proteins: a 5' uridine containing microRNA loads into AGO1, while one with a 5' adenosine loads into AGO2, and one with cytidine at the 5' terminus loads into AGO5 [6]. Furthermore, the strong 5' bias for uridine also suggests selective loading of the mature strand microRNA over the star strand. This is consistent with the hypothesis that only one of the strands enters into the RISC complex [6]. Our datasets were generated from total RNA that was isolated from the whole body of adult H. vitripennis and hence may lack the specificity, depth, and resolution that is needed for more robust analysis of temporal and spatial expression of the microRNAs. Since the focus of our study was the adult stage of H. vitripennis, the expression pattern of the microRNAs in the juvenile and embryonic stages was not determined. Conserved microRNAs are dominant in late stage embryos of D. melanogaster and D. virilis, whereas in the early embryonic stages rapidly evolving microRNAs are more common [52]. Thus, both the conservation (i.e., whether conserved or novel) and type (i.e., whether mature or star strand) of microRNA that is present at a particular time can differ depending upon the developmental stage of the organism. Our selection of only a single developmental stage of H. vitripennis may partly explain the generally low occurrence of star strand reads in our library (Table 5). Recent studies also suggest that the presence of star strands is more common in eukaryotes and vertebrate organisms [5] and that these star strands may function as mature microRNAs [53]. Consistent with these studies is the report of the enhanced abundance of certain star strands in microRNA libraries from Manduca sexta [30] and the zebra finch [53]. It will be interesting to characterize the novel microRNAs of H. vitripennis with respect to the developmental stage and under different physiological conditions. Interestingly, the putative origin of 7 out of the 14 identified novel microRNAs (Hvi-miR13059, Hvi-miR9237, Hvi-miR19117, Hvi-miR24402, Hvi-miR41359, Hvi-miR91 and Hvi-miR24657) was found within a transposable element sequences. Insect genomes are substantially loaded with transposable elements, and their presence helps to explain the varying sizes of their genomes. In particular, insect genomes carry widely varying amounts of repeat sequences in their genomes [54]. In H. vitripennis, we previously identified four different types of transposable elements with a majority of them belonging to TC3 type (131 copies) and piggy-Bac elements (58 copies) [15]. Recent evidence suggests that the functional evolution of micro-RNAs is partly due to insertions by the transposable elements resulting in the generation of novel microRNAs [55][56][57].
MicroRNAs influence the development and immunity of the host organism by regulating mRNA levels [48][49][50][51]. The majority of microRNAs in animals are shown to interact with their target mRNAs through the 3'-and 5'-untranslated regions (UTRs) of the mRNAs. MicroRNAs are thus believed to regulate their target genes at a post-transcriptional level. We used a library of full-length transcripts (~52,700 putative mRNAs) of H. vitripennis for in silico target prediction.  We identified 9,132 potential microRNA targets through the Miranda target prediction program (Tables 6 and 7, S2 and S3 Tables). Several key regulatory genes, as well as genes involved with signaling pathways, ion channels, enzymatic processes, catabolic reactions and KEGG pathways were identified by this process. The future validation and further study of these potential targets will lead to a deeper understanding the regulatory biology of H. vitripennis and likely lead to the identification of novel gene targets that can be exploited for control of this pest insect.

Conclusions
In the present study, we identified 14 novel microRNA candidates from H. vitripennis. To date, this is the first report on such an analyses in H. vitripennis. These data provide a promising approach in gene functional studies through the use of RNAi-based approaches in cell lines or whole insects. The expression of microRNAs that target key regulatory genes in H. vitripennis through transgenic plants has potential to control H. vitripennis and/or the pathogens that it vectors.
Supporting Information S1