Identification of Aedes aegypti Long Intergenic Non-coding RNAs and Their Association with Wolbachia and Dengue Virus Infection

Long intergenic non-coding RNAs (lincRNAs) are appearing as an important class of regulatory RNAs with a variety of biological functions. The aim of this study was to identify the lincRNA profile in the dengue vector Aedes aegypti and evaluate their potential role in host-pathogen interaction. The majority of previous RNA-Seq transcriptome studies in Ae. aegypti have focused on the expression pattern of annotated protein coding genes under different biological conditions. Here, we used 35 publically available RNA-Seq datasets with relatively high depth to screen the Ae. aegypti genome for lincRNA discovery. This led to the identification of 3,482 putative lincRNAs. These lincRNA genes displayed a slightly lower GC content and shorter transcript lengths compared to protein-encoding genes. Ae. aegypti lincRNAs also demonstrate low evolutionary sequence conservation even among closely related species such as Culex quinquefasciatus and Anopheles gambiae. We examined their expression in dengue virus serotype 2 (DENV-2) and Wolbachia infected and non-infected adult mosquitoes and Aa20 cells. The results revealed that DENV-2 infection increased the abundance of a number of host lincRNAs, from which some suppress viral replication in mosquito cells. RNAi-mediated silencing of lincRNA_1317 led to enhancement in viral replication, which possibly indicates its potential involvement in the host anti-viral defense. A number of lincRNAs were also differentially expressed in Wolbachia-infected mosquitoes. The results will facilitate future studies to unravel the function of lncRNAs in insects and may prove to be beneficial in developing new ways to control vectors or inhibit replication of viruses in them.


Introduction
Dengue and Zika viruses are related mosquito-borne viruses that have a common primary vector, Aedes aegypti and infect millions of people worldwide [1,2]. Recent outbreaks of Dengue and Zika in South America pose a serious risk for other tropical regions in the world as Ae. aegypti is one of the most abundant mosquito species in these areas [2]. Although certain vaccines have been licensed in some countries, there are no efficient specific therapeutics available for either diseases, hence, the best protection against their global spreading is an efficient vector control program [3,4].
The genome sequence of Ae. aegypti is available, however, it has not been fully annotated. Only 2% of its large genome (1.376 Mb) has been annotated as protein coding genes and it reflects the presence of great proportions of non-coding transcripts as well as repetitive elements [5]. Transcriptomic changes, including those of non-coding transcripts, could provide a genome scale insight into host-pathogen interactions. Previous studies identified a series of small ncRNAs in Ae. aegypti and demonstrated their interaction with arboviruses [6][7][8][9], but our knowledge about their long ncRNAs is limited.
RNA transcripts longer than 200 nucleotides, which do not contain an open reading frame of longer than 100 amino acids, are simply defined as long ncRNA [10]. Generally, they are classified by their location relative to their neighboring protein-coding genes and include the long intergenic ncRNA (lincRNA), intronic lncRNA, antisense lncRNA and enhancer RNA [10]. Although a number of mammalian lncRNAs have been characterized and identified in the last few years, genome-wide identification of this class of ncRNAs has only recently become possible with the arrival of deep sequencing technologies. An expanding body of evidence reveals that lncRNAs, once described as dark matter, are involved in many biological processes such as genomic imprinting and cell differentiation [11]. They also play important roles in epigenetic and non-epigenetic based gene regulation [12]. Relatively, little is known about their involvement in activation and differentiation of immune cells, but new discoveries have revealed the involvement of lncRNA in defense systems [13]. Previous works have also outlined their quick responses to different stimuli and stress factors [14][15][16][17]. In addition, it has been shown that some lncRNAs enhance virus replication or decrease antiviral immunity [18].
Although in most host-virus interaction studies typically protein-coding genes have been the center of attention, there are few examples of virus and host lncRNA interactions in human and mouse models [18,19]. For instance, Hepatitis B virus (HBV) infection altered lncRNA profiles in patients, with about 4% of human lncRNAs showing more than 2-fold changes in HBV infected liver tissue [20]. Winterling et al. (2014) identified a virus inducible lncRNA, which is induced by vesicular stomatitis virus and several strains of Influenza A virus (IAV) [18].
The sequence and structure of lncRNAs are important in their function, in particular for their interaction with DNA, RNA, or proteins. In case of extensive base-pairing of lncRNA with target mRNA, translation can be stabilized, while partial base-pairing may accelerate mRNA decay or inhibit translation of the target mRNA [21]. It has been shown that some lncRNAs interact with other small ncRNAs such as miRNAs. For example, in silkworm, 69 lncRNAs originating from 33 gene loci, may serve as miRNA precursors, and 104 lncRNAs may function as competing endogenous RNAs (ceRNAs) [22]. LncRNAs are also targeted by miRNAs similar to mRNAs and reduce their stability. They may also act as sponge or decoy of miRNAs, and compete with miRNAs for binding to mutual target mRNAs [21].
In insects, only a few genes have been experimentally annotated as lncRNA. It has been estimated that more than 5000 loci potentially encode non-coding transcripts in Drosophila melanogaster, however, just seven loci (bxd, Hsrω, pgc, roX1, rox2, sphinx and yar) have been annotated as functional regulatory lncRNAs by experimentally derived data [23,24]. We recently found that a number of lncRNAs in Plutella xylostella, a pest of cruciferous plants, were linked to the insect's resistance to insecticides and might be involved in detoxification processes [14]. Jenkins et al (2015) identified 2,949 lncRNAs in the malaria mosquito vector, Anopheles gambiae, using RNA-Seq data [25]. They showed that in various Anopheles species, lncRNAs have considerably lower sequence conservation as compared with protein-coding genes. In another study, it has been shown that 43% of total midgut transcripts of An. gambiae are lncRNAs and 32% of them showed some level of homology to other species [26].
The current study generated a comprehensive list of Ae. aegypti lincRNAs, which will be a complement to the other ncRNAs (microRNAs and piRNAs) that have already been discovered in this medically important species. This work also helps to improve the present annotation of the genome of Ae. aegypti. We also examined the expression pattern of some selected lincRNAs in response to microbial challenge namely dengue virus serotype 2 (DENV-2) and Wolbachia infections to identify potential immune-related lincRNAs in Ae. aegypti [27,28]. The results help better understanding of mosquito-pathogen interactions providing new insights on the potential role of lncRNAs as candidates for exploitation to inhibit replication of mosquito-borne viruses.

RNA-Seq Data preparation
Previously sequenced RNA-Seq raw data of Ae. aegypti were downloaded from NCBI Sequences Read Archive and ArrayExpress Archive with accession numbers SRA048559, SRA058076, SRA244067 and ERP002530 [29][30][31][32]. Raw data were stripped of adapters using CLC Genomic Workbench version 7.5.1 and reads with quality score of above 0.05 and maximum 2 ambiguous sequences were retained for further analysis.

Large gap mapping and transcript discovery
The CLC Genomic workbench's Transcript Discovery plugin was used for lincRNA discovery in the Ae. aegypti genome. New transcripts were identified by large gap mapping of 1,148,814,115 reads of 35 RNA-Seq libraries to the genomic reference (AaegL3.3). We implemented strict mapping criteria (mismatch, insertion and deletion costs: 2: 3: 3 respectively). The minimum similarity and length fraction of 0.9 between a mapped segment and the reference were allowed as part of the mapping criteria. The large gap mapper algorithm also requires each mapped segment to include at least 10% of the read with minimum length of 17 bases. We considered a gap with maximum of 50 Kbp distance between mapped read segments to span the introns from RNA-Seq data. The annotations were generated by inspecting mapping of reads and identifying likely regions corresponding to genes, including their exons and splice sites. The algorithm scans each gap in the read mapping to explore whether the gap is assigned to a valid splice site or can be relocated to a valid splice site without cost.

lincRNA identification pipeline
A rigorous filtering pipeline was developed to remove transcripts that may potentially encode proteins. The pipeline for Ae. aegypti lincRNA discovery is summarized in Fig 1. We identified 75,069 potential genes using the CLC Genomic Workbench transcript discovery algorithm. The genes that were annotated as known Ae. aegypti protein-coding genes were discarded and 30,865 potential genes were also checked for any exon or intron overlap with other known Ae. aegypti protein-coding genes. We selected 22,079 sequences, which were located more than 1kb away from any other known transcripts, for finding putative open reading frames (ORF). All possible six frames were produced for all selected sequences and then the translated sequences were subjected to a domain search to identify any putative conserved protein domains through Pfam v27.0 database [33]. We discarded 8,795 sequences with potential ORF above 100 aa or conserved protein domains. The remaining sequences were submitted to a coding potential assessment tool (CPAT), which utilizes a logistic regression model made with four sequence features: ORF size, ORF coverage, Fickett TESTCODE statistics and hexamer usage bias [34]. We applied the coding probability threshold of 0.3, which led to discarding 376 sequences as putative coding RNAs. We also implemented an expression threshold on our data to strengthen the identification pipeline. Sequences with more than 10 mappable reads in at least 17 out of 35 RNA-Seq libraries were considered as valid sequences and were kept for the next step. Any possible similarity with other known proteins was found by using BLASTx algorithm against nr and Swiss port database (E-value cut off 10 −5 ). Finally, 3,842 potential lincRNAs were identified and used for further analysis.
To identify Ae. aegypti putative lincRNAs that are regarded as small RNA associated lincR-NAs, we used the Blast algorithm to search for Ae. aegypti precursor miRNA sequences in the predicted Ae. aegypti lincRNA dataset. We also used publicly available small RNA libraries from DENV-infected and non-infected samples (SRP026241) in this analysis for further characterization of lincRNA_1317. All known Ae. aegypti miRNA sequences were mapped to lincRNA_1317 for possible best fitting using RNAhybrid, which is a tool for finding the normalized minimum free energy (mfe) of RNA. We did not allow G:U pairing in the seed region (nucleotide 2-8) and required miRNA-lincRNA duplexes to have a helix in this region. Maximum 5nt were allowed as unpaired nucleotides in either side of an internal loop. LncTar algorithm [35] was used to explore any potential interaction between lincRNA_1317 and DENV-2 genome (accession no. NC_001474.2) by finding the normalized mfe joint structure of two RNA molecules based on base pairing.

Identification of differentially expressed lincRNAs upon infection
The Ae. aegypti genome was annotated with the final list of lincRNAs and used as reference for RNA-Seq analysis in CLC Genomic Workbench. To measure the lincRNA normalized expression value, RPKM (Reads Per Kilobase per Million reads) was assigned for each library [36]. To find the differential expression pattern of lincRNAs in response to DENV infection, data from DENV-2 (Jam1409) infected midgut and carcass tissues at 4 days post-infection (dpi) were compared with their corresponding control groups [30]. Baggerley's test, a count based statistical analysis was done on the data. The samples were given weights depending on their total counts. Based on the test the weights are obtained by supposing a Beta distribution on the proportions in a group, and estimating these, along with the proportion of a binomial distribution. We selected 20 potential lincRNAs with more than 4-fold change for further analysis with RT- qPCR in DENV-2 (New Guinea C strain) infected Ae. aegypti cell line (Aa20) and screened their expression profile in Wolbachia-infected mosquitoes.
Ae. aegypti infected with the wMelPop-CLA strain of Wolbachia (+Wol) and without Wolbachia (-Wol, tetracycline-cured line) were stocks produced previously [37]. For the experiments in this work, 4-day-old female mosquitoes were used from which total RNA was extracted with 6-10 adult mosquitoes for each biological replicates.
Expression analysis of Ae. aegypti lincRNAs Detection and validation of the relative abundance of selected lincRNAs was carried through lincRNAs' specific primers using SYBR Green chemistry in real time PCR machine. Briefly, total RNA was extracted from cells using Qiazol reagent according to the manufacturer's instructions (Qiagen). The TURBO DNA-free kit (Ambion, USA) was used to remove possible genomic DNA contamination in RNA samples. First strand cDNA was synthesized from 2.5 μg of RNA using a poly-dT primer and Superscript III reverse transcriptase (Life Technologies). qPCR primers were designed using primer design tool of NCBI [38]. QuantiFast SYBR Green PCR Master Mix with ROX was used to quantify the relative expression of lincRNAs between different treatments. Three independent biological replicates were considered along with three technical replicates for each treatment. Reactions were performed in a Rotor-Gene thermal cycler (Qiagen) under the following conditions: 95°C for 5 min, and 40 cycles of 95°C for 10s and 60°C for 30s, followed by the melting curve (68°C to 95°C). Melting curves were analysed to examine the specificity of amplification. Relative expressions were calculated using the Rotor-Gene software and the mosquito RPS17 as reference gene for normalization. Unpaired t-test was used to identify statistically significant differences.

RNAi of selected lincRNAs and virus replication assay
To check the functional importance of the identified novel lincRNAs, dsRNAs were synthesized to knockdown selected lincRNAs (2329, 1613 and 1317) to check their effect on DENV replication. Briefly, primers with added T7 promoter sequence (S1 Table) were used to generate 250-600 bp PCR products from selected lincRNAs. Megascript T7 kit (Ambion) was used according to the manufacturer's instruction to generate respective dsRNAs. To induce efficient RNA silencing, Ae. aegypti Aa20 cells were double transfected with dsRNAs against selected lincRNAs. Aa20 cells were re-suspended and~3×10 5 cells were added to each well of a 12-well plate. Cells were allowed to settle for~1 h, medium was removed and replaced with a transfection mixture consisting of 0.5 ml medium (1:1 Schneider medium and Mitsuhashi-Maramorosch with 10% FBS), 8 μl Cellfectin (Invitrogen), and 5 μg dsRNA either for selected lincRNAs or GFP as control. Cells were also treated with 3 μg dsRNA 72 h after the primary transfection to increase the silencing efficiency of selected lincRNAs. Six hours after the secondary transfection, cells were infected at 1 multiplicity of infection (MOI) with DENV-2 (New Guinea C strain). All the treatments were collected three days post-infection. RNA extraction and cDNA synthesis were carried out as above. qPCR was performed to confirm the knockdown and the effect of particular lincRNA knockdown on the genomic RNA of DENV-2. Each treatment was repeated three times. All data from three biological replicates were subjected to one-way ANOVA statistical analysis. Brown-Forsythe test was used to check the equality of group variances and Tukey's multiple comparisons test was also used to examine significant statistical differences among treatments.

Results and Discussion
Identification and characterisation of Ae. aegypti lincRNAs In total, 3,482 putative lincRNAs in 1,114 Ae. aegypti genome scaffolds were identified (S2 Table). The Ae. aegypti lincRNA genes displayed a slightly lower GC content (mean: 40.1%) in comparison to 47.8% in their protein-coding gene sequences (Fig 2A). The lower GC content or AT enrichment is a typical characteristic of lincRNAs and our findings are congruent with predicted lincRNAs in other species [14,39,40]. The majority of Ae. aegypti predicted lincRNAs are smaller than 3000 bases and their length distribution is represented in Fig 2B. These mosquito lincRNA candidates are notably shorter in length than protein-coding genes, demonstrating another well-known characteristic of lincRNA transcripts (Fig 2C) [41,42]. The majority of Ae. aegypti genome scaffolds contain less than five lincRNA loci (~80%), however, 23 of scaffolds (2%) were enriched with more than 10 lincRNAs (Fig 2D). The detailed information of these scaffolds, which contain the highest number of lincRNAs are summarized in Table 1.
We examined all the identified lincRNAs to determine their association with Ae. aegypti miRNA precursors and piRNA clusters. We found that the pre-miRNA sequences of aae-miR-2940 and aae-miR-285 are located in lincRNAs 1431 and 3299, respectively. We could not detect any other pre-miRNA sequences identified in Ae. aegypti in the lincRNAs. Also, lincRNA 1978 and 792 are originated from two previously reported piRNA clusters [43] located at supercontig 1.478 and 1.98, respectively.
LincRNAs demonstrate low evolutionary sequence conservation even among closely related species [10,14]. We used the BLAST algorithm bit score to identify the level of similarity among Ae. aegypti lincRNA sequences with other closely related insect genomes such as Aedes albopictus, Culex quinquefasciatus and Anopheles gambiae ( Fig 3A). As expected, most of the identified lincRNAs showed high level of similarity with Ae. albopictus genome sequence and probably are genus specific. The E-value cut off 10 −50 was applied to our screening with the BLAST algorithm to identify the conserved sequences. Although the Ae. aegypti lincRNAs shared high level of sequence similarity with the genome of Ae. albopictus, only 62 and 7 lincR-NAs had sequence similarity with Cx. quinquefasciatus and An. gambiae, respectively (Fig 3B). They were mostly limited to a single short region with high conservation.

Ae. aegypti lincRNAs change upon microbial challenge
Following the identification of Ae. aegypti lincRNAs, we analyzed their transcript levels in DENV-2 infected mosquitoes. To produce the lincRNA profile of infected and non-infected mosquitoes, we re-analyzed previously published RNA-Seq data from Ae. aegypti midgut and carcass samples at 4 dpi [30] (Fig 4). 248 and 203 lincRNAs with fold changes above four were identified in the RNA-Seq libraries of midgut and carcass, respectively (S3 Table). The majority of differentially expressed lincRNAs were considerably overexpressed in both samples. The abundance of only 32% of Ae. aegypti lincRNA candidates decreased in response to DENV-2 infection in the mosquito carcass sample. Thirty lincRNAs were differentially expressed in both examined samples. The transcription levels of 72 lincRNAs increased after infection while their expression could not be detected in the non-infected midgut tissue sample. We selected 20 candidates of those differentially expressed lincRNAs from RNA-Seq analysis data for further investigation. The relative expression of lincRNA candidates were examined by reverse transcription quantitative polymerase chain reaction (RT-qPCR) upon DENV-2 infection in Aa20 mosquito cells. Only significantly overexpressed lincRNAs after DENV-2 infection are represented in Fig 5. Although we used Aa20 cells for the lincRNA expression assays, the expression patterns of almost all the examined lincRNAs (5 out of 6) were consistent with the RNA-Seq data (adult mosquito carcass sample). We used a poly-dT primer to produce cDNA, which also confirmed that all of those identified transcripts have poly-A tails and therefore are true transcripts. Based on these results, significant increase in the transcription levels of a selected number of Ae. aegypti lincRNAs suggests their possible involvement in host-pathogen interaction but further investigations are required to confirm their roles in antiviral/ immune responses.
We also examined the impact of an endosymbiotic bacterium, Wolbachia, on some selected Ae. aegypti lincRNAs, which showed significant changes in response to DENV-2 infection. This gram-negative bacterium is transmitted maternally and potentially infects more than 40% of all insect species, manipulating its hosts using different strategies [44,45]. A fascinating aspect of Wolbachia infection is limiting replication of vector-borne pathogens in mosquitoes [45,46]. However, the mechanism(s) behind virus blocking is largely unknown. Here, we found that the transcript levels of several lincRNA genes significantly increased in Wolbachia-infected Ae. aegypti mosquitoes (Fig 6), which may lead to differential regulation of cellular proteincoding genes. Our previous studies showed that Wolbachia could manipulate host small ncRNAs such as miRNAs and piRNAs [47]. An overall induction of small ncRNAs between 18 and 28 nucleotides was also observed in Ae. aegypti cell line infected with wMelPop-CLA strain of Wolbachia [48]. It was assumed that the upregulation of small ncRNAs in infected cells may result in an enhanced immune response and activated RNAi pathway. However, the role of these modifications in the host lincRNA gene expression profile, and potentially in anti-viral responses, is unknown and may lead to the discovery of lincRNAs that could be utilized for inhibition of virus replication in mosquitoes.
A recent study on mouse bone marrow-derived macrophage (BMDM) model reported a significant upregulation in 72 lincRNAs after treatment with the synthetic bacterial lipoprotein Pam3CSK4, which acts through Toll-like receptor [49]. In another study, differential expression of approximately 500 annotated mouse lncRNAs was reported during infection with severe acute respiratory syndrome coronavirus [50]. Recently, it has been shown that honeybee lincRNAs are also differentially expressed during infection with various viruses such as sacbrood virus (SBV) and deformed wing virus (DWV), but the biological significance of these lincRNAs is completely unknown [51]. Although exploring the in vivo functions of immunerelated lincRNAs is one exciting area for future studies, the differential expression of some lincRNAs could simply be byproducts of mRNA biogenesis or changes in global transcriptional profile due to microbial challenges [52,53]. Struhl (2007) believed that the transcriptional machinery is not perfect producing RNAs that serve no purpose or have no significant role in infection [54]. On the other hand, there are several examples which have shown that lincRNAs could be potentially important factors in host antimicrobial responses, and may represent a new class of signaling molecules involved in innate immunity or provide a new layer in gene regulation. For instance, two interferon (INF) induced lncRNAs, which were upregulated by influenza and vesicular stomatitis viruses, regulate the expression of the antiviral factor tetherin in human HuH7 cells [55].

RNAi of selected lincRNAs and enhancement of DENV-2 replication
To confirm the role of DENV-induced lincRNAs in viral replication, we used RNAi-mediated silencing of two selected lincRNAs (lincRNA_1317 and 1613) using dsRNA in Aa20 cells followed by DENV-2 infection. Only RNAi-mediated silencing of lincRNA_1317 led to enhancement of DENV-2 replication (Fig 7A). Silencing of the lincRNA was confirmed by RT-qPCR ( Fig 7B). Interestingly, expression of Ae. aegypti lincRNA_1317 increased substantially following the progression of infection (Fig 7C) suggesting that this lincRNA might be involved in antiviral response. This idea is consistent with the finding that lincRNA_1317 was also highly overexpressed (2.33 fold) in Wolbachia-infected mosquitoes as compared with non-infected mosquitoes (Fig 6).
While there are no reports on the involvement of lncRNAs in host-pathogen interactions in insects, time-dependent over-expression of host lincRNAs in response to viral infection has been observed in humans. A recent study showed more than 80% of host cell lncRNAs were upregulated upon an adenovirus infection of human primary lung fibroblast cells [56]. Zhang et al. (2013) reported alterations of expression of cellular lncRNAs in HIV-1-infected T cells. Among differentially expressed lncRNAs, NEAT1 expression notably increased in infected cells. When NEAT1 was silenced, virus production was enhanced by increasing the nucleus-tocytoplasm export of HIV-1 transcripts containing Rev-dependent instability element [57]. A significant induction in this lncRNA expression in response to influenza virus and herpes simplex virus infection has also been shown [58].
To further investigate the potential role of lincRNA_1317 in mosquito-pathogen interaction, we determined its association with host endogenous small RNAs and its possible direct interaction with DENV. Although this lincRNA is not located in any of the known piRNA clusters, the majority of mappable small RNA reads to its sequence are in the range of 26-29 nt (S1 Fig). However, there was no difference in the mapping pattern and mapped read length distribution when reads from DENV-infected and non-infected small RNA libraries were mapped to lincRNA_1317 (S1 Fig). It has been shown that piRNA-like small RNAs have a large impact on lincRNA transcriptome [57], but our knowledge about the function of piRNA-mediated lncRNAs is still limited. Recently, it has been reported that piRNAs derived from transposons and pseudogenes facilitate the degradation of lncRNAs in mouse late spermatocytes [57].
Next, we hypothesized that Ae. aegypti lincRNA_1317 response to microbial challenge could be due to cross-regulation between miRNAs and the lincRNA. Ae. aegypti miRNA recognition elements on lincRNA_1317 were identified by calculating the normalized minimum free energy (mfe) of hybridization for each Ae. aegypti miRNA and lincRNA_1317 using RNAhybrid core script. Binding site enrichment was detected for a few miRNAs with more than two recognition elements (Table 2). For instance, more than four recognition sites were predicted for miR-278-5p and miR-252-3p on lincRNA_1317. We also identified some hot spots for miRNA recognition sites on lincRNA_1317, which may allow multiple miRNAs to bind to the same regions (S2 Fig). miRNAs can reduce lincRNA stability by targeting their transcripts similar to mRNAs. Also, lincRNAs with multiple recognition sites may actually be competitive inhibitors of miRNA function and stopping them from binding to their genuine targets by sequestering them [21]. Although the mfe for some of those miRNA-lincRNA recognition sites suggests high probability of a binding event, further experimental investigations are required to validate this interface.
We also used LncTar algorithm to predict any direct interaction between lincRNA_1317 and DENV-2 genome. One potential interaction was predicted in the region 1-3370 of lincRNA_1317 and the region of 3210-6579 of DENV-2 genome with mfe of -61.73 (normalized dG -0.0184). This tool has accuracy rate of 80% [35], but does not consider the tertiary Possible involvement of lincRNA_1317 in DENV-2 replication in Ae. aegypti Aa20 cells. A) Aa20 cells were double transfected with the transfection reagent only (Mock), dsRNA to GFP (dsCont) or dsRNA to lincRNA_1317 (dslncRNA) for three days followed by 1 MOI infection of the cells with DENV-2. Primers to the NS2A region (S1 Table) were used for measuring the relative DENV-2 genomic RNA levels. B) RNAi silencing of lincRNA_1317 using dsRNA was confirmed by RT-qPCR. C) Changes in DENV-2 genomic RNA levels during the course of infection analyzed by RT-qPCR on RNA extracted from Aa20 cells at 1, 3 and 5 days post-infection (dpi). D) DENV-2 infection increased the transcript levels of Ae. aegypti lincRNA_1317 in 3 rd and 5 th days post-infection when RNA from cells were analyzed by RT-qPCR. In all the experiments shown in this figure, three biological replicates, each with three technical replicates were used. **, p < 0.01; ***, p < 0.001. doi:10.1371/journal.pntd.0005069.g007 structure of RNA, which could play a role in RNA-RNA interactions, and therefore further studies are required to validate any potential interaction.
The involvement of lincRNA_1317 in host response to viral infection might be through its interactions with regulatory proteins that are involved in epigenetic changes by directly interacting with chromatin modifying enzymes or DNA binding proteins such as transcription factors. This interaction has been shown in several examples in mammalian systems, including host-virus interactions in which lncRNAs mediate antiviral responses by controlling the expression of immune-related genes (reviewed in [58]).
Although our knowledge of the biological function of this class of ncRNAs in mosquitoes is still limited, the results generated from this study will facilitate forthcoming explorations of lincRNA functions in insects. Clearly, further research is required to provide concrete experimental evidence to support the role of lincRNA_1317 or any other Ae. aegypti lincRNAs in host-pathogen interaction. With advances in technology, the mosquito lincRNA-protein interactions can be identified using high-throughput sequencing of immunoprecipitated RNA after cross-linking (CLIP-Seq). Further, functional studies could be carried out to characterize immune-related lincRNAs. The involvement of lincRNAs in pathways associated with responses to viral infection and cellular stress makes them interesting candidates as potential targets for manipulation to inhibit virus replication or control vector populations.