Expression Profiling Analysis Reveals Key MicroRNA–mRNA Interactions in Early Retinal Degeneration in Retinitis Pigmentosa

Abstract

Purpose: The aim of this study was to identify differentially expressed microRNAs (miRNAs) that might play an important role in the etiology of retinal degeneration in a genetic mouse model of retinitis pigmentosa (rd10 mice) at initial stages of the disease.

Methods: miRNAs–mRNA interaction networks were generated for analysis of biological pathways involved in retinal degeneration.

Results: Of more than 1900 miRNAs analyzed, we selected 19 miRNAs on the basis of (1) a significant differential expression in rd10 retinas compared with control samples and (2) an inverse expression relationship with predicted mRNA targets involved in biological pathways relevant to retinal biology and/or degeneration. Seven of the selected miRNAs have been associated with retinal dystrophies, whereas, to our knowledge, nine have not been previously linked to any disease.

Conclusions: This study contributes to our understanding of the etiology and progression of retinal degeneration.

Retinitis pigmentosa (RP) is a heterogeneous collection of inherited retinal degeneration diseases that leads to vision impairment and for which there is no standard treatment. It is characterized by the progressive death of retinal photoreceptor cells. The primary defect underlying RP affects the function of photoreceptors and/or retinal pigment epithelium cells, where molecular and cellular mechanisms trigger the apoptotic degeneration of rods and in many cases cones. 

Despite there being more than 350 genes and/or loci associated with retinal dystrophies, the precise mechanisms that link mutations in these genes with retinal degeneration are still not well defined. In addition to many molecular pathways involved in both survival and death of photoreceptors and RPE cells, there is growing evidence of the role of epigenetic mechanisms as key players in gene expression regulation, retinal development and function, and photoreceptors survival.^1,2 These include effects mediated by microRNAs (miRNAs), small RNAs that partially inhibit the translation of a transcript or cause mRNA degradation.³ According to the current version of the miRBase database (www.mirbase.org), more than 2000 human miRNAs have been identified, each of which targets several hundreds of different mRNAs, mostly due to the fact that most target sites on the mRNA have only partial base complementarity with their corresponding miRNA. Moreover, individual mRNAs may contain multiple binding sites for different miRNAs, resulting in a complex regulatory network. Taken together, it appears that relatively few miRNAs can regulate as many as 50% to 60% of genes in humans and other mammals.⁴ 

Control of protein expression mediated by miRNA modulation has been described as a widely used mechanism for posttranscriptional regulation of metabolic pathways.^5,6 There is growing evidence of the involvement of miRNAs dysregulation in a broad spectrum of health problems, including various types of cancer and heart and inflammatory diseases, as well as in hereditary diseases such as cystic fibrosis.^7–15 

Current molecular tools have allowed the selective modulation of miRNA expression. This has led to the development of what will become the first clinical trials based on an antisense approach seeking to regulate miRNA function in patients with cardiovascular conditions.¹⁶ 

More than 250 miRNAs have been identified in the retina, and recent studies have shown the relevance of Dicer, a key enzyme in miRNA biogenesis, in photoreceptor survival.^2,17 Similar to many other human conditions, retinal diseases also involve miRNA dysregulation. Indeed, there is a growing list of miRNAs that are known to be differentially expressed in retinal disorders in humans and in animal models, including miR-9, miR-34a, miR-125b, and miR-155 in macular degeneration,^18,19 miR-146a and miR-195 in diabetic retinopathy,^20,21 and miR-125a and miR-17 in retinoblastoma.^22,23 In addition, a group of antiapoptotic miRNAs (miR-155, miR-146a, and miR-29b) was found to be differentially expressed during the chronic cell death stage in canine models of retinal degeneration.²⁴ Furthermore, a common pattern of aberrantly expressed miRNAs was described in four murine models of RP where miR-1, miR-133, and miR-142 were upregulated despite different genes (Rho and Rds) and inheritance patterns being involved.^25,26 

To identify potential miRNAs associated with the etiology of RP, we used the retinal degeneration 10 (rd10) mouse: a well-established mouse model of autosomal recessive RP, which carries a spontaneous point mutation in exon 13 of the β subunit of the rod phosphodiesterase gene (Pde6β). This mutation results in deregulation of intracellular calcium homeostasis that leads to photoreceptor cell death via apoptosis.^27–29 Mutations in the PDE6β gene also cause RP in humans accounting for approximately 5% of all autosomal dominant RP cases worldwide.^30,31 Therefore, rd10 mice provide a useful and widely used model of RP that closely resembles the human disease emulating the progression of typical autosomal recessive RP. 

In summary, we analyzed the expression pattern of miRNAs and associated mRNAs during the early stages of retinal degeneration. Our results provide valuable insights into the molecular mechanisms underlying RP, including the identification of differentially expressed miRNAs, and the association of predicted mRNA targets in the retina of the rd10 mice. These include both miRNAs previously reported in retinal studies, such as miR-1a, miR-133a, miR-146a, miR-17, miR-29b, and miR-155, and miRNAs not previously ascribed to any disease, namely miR-6937, miR-6240, miR-3473b, and miR-6970. This study might also help to identify potential miR-based therapeutic approaches. 

Materials and Methods

Ethics Statement and Animal Handling

Animal handling and experiments were conducted in accordance with the ARVO Statement for the Use of Animals in Ophthalmic and Vision Research and were approved by the Animal Care and Use Committee of Donostia University Hospital and the Clinical Research Ethics Committee of the Basque Country, Spain. Studies were performed using age-matched wild-type (WT) C57BL/6J mice as controls and a congenic inbred strain of rd10 mice, both purchased from the Jackson Laboratory (Bar Harbor, ME, USA). Animals were maintained under a 12-hour light/12-hour dark cycle at 22°C, with controlled humidity, ranging between 45% and 55%, and with water and food provided ad libitum at the Animal Facility of BioDonostia Health Research Institute in San Sebastian (Spain). 

Sample Collection

Rd10 and age-matched WT mice were euthanized by cervical dislocation following CO₂ inhalation. All mice were euthanized at the same time of day to minimize possible bias from circadian effects on miRNA expression. 

Retinal tissue from each animal was isolated and kept at 4°C within a few minutes until use. Total RNA isolation was performed using the miRNeasy Mini Kit and RNase-Free DNase Set (Qiagen, Woburn, MA, USA), in accordance with the manufacturer's protocol. RNA was stored at −80°C until use. Only those RNA samples with A260/A280 ratios greater than 1.8 were used, as determined by a Nano-Drop ND-1000 spectrophotometer (Nanodrop Technologies, Wilmigton, DE, USA). Furthermore, sample quality was tested in each of the different RNA expression profiling arrays described below. 

Eyes were enucleated, fixed in PBS containing 4% paraformaldehyde, and cryoprotected in PBS containing 30% sucrose. They were then embedded in OCT (Tissue Tek, Sakura Finetek, Tokyo, Japan) and sectioned at 7 μm in a cryostat. 

Time Course of Photoreceptor Cell Death

We used a terminal deoxynucleotidyl transferase-mediated dUTP nick-end labeling (TUNEL) assay to study the photoreceptor cell death process to select the time points for collection of retinal samples for our miRNA/mRNA expression studies. We analyzed retinal sections from postnatal days 13 to 22 (P13 to P22), obtained at 12-hour intervals using the DeadEnd Fluorometric TUNEL System (Promega, Madison, WI, USA). At least three samples were analyzed for each experimental group and time point. 

Global miRNA and mRNA Expression Profiling

A detailed global analysis of miRNA expression and its correlation with the transcriptome was performed using the miScript miRNA PCR array (Qiagen), GeneChip miRNA 4.0 array and GeneChip Mouse Transcriptome Array 1.0 (MTA 1.0) (Affymetrix, Santa Clara, CA, USA), covering 100% of miRBase v20 mouse miRNAs and 100% of the mouse transcriptome. 

Three miScript miRNA PCR Arrays (Qiagen; MIMM-001Z) were used to analyze 12 samples from pooled retinas, corresponding to two samples per experimental group, rd10 and WT, from three time points: P13, P15, and P17. Each sample consisted of a pool of nine retinas from six mice. 

Ten GeneChip miRNA 4.0 arrays and 10 GeneChip Mouse Transcriptome Arrays 1.0 were used to analyze 10 samples from pooled retinas: 3 for each rd10, from P13, P15, and P17 and 1 for WT from P13. Each sample consisted of a pool of nine retinas from six mice (Supplementary Fig. S1). 

We first performed analysis of raw data from all GeneChip arrays, including a detection step (a probe set being detected above background with an associated P value) resulting in an absence/presence call and a quantile normalization step using the Expression Console Software v 1.4.1 (Affymetrix). In a subsequent filtering step, all nonmouse probe sets and those that were not detected (values below the background signal) in any of the arrays were removed to ensure that we worked with miRNAs and mRNAs actually expressed in mouse retina. 

For identification of differentially expressed (DE) miRNAs and mRNAs, we used as reference rd10 and WT mice at P13 (in GeneChip arrays), 3 days prior to the onset of apoptosis. Adjusted P < 0.05 was considered statistically significant, and fold-change cutoffs of ±1.5 and ±1.3 were used for miRNAs and mRNAs, respectively. Further criteria used included the exclusion of miRNAs with values close to the background signal and triplicates with a high SD. 

We then filtered the list of DE miRNAs based on fold-change values and on initial prediction studies of their putative target mRNAs (i.e., whether the corresponding genes were potentially relevant to retinal pathways in health and disease, such as genes known to be involved in apoptosis, inflammation, or normal retinal function or previously reported to be related to RP). 

Candidate miRNAs were further classified into two groups based on their expression dynamics (i.e., whether altered expression occurred prior to or just after the onset of photoreceptor apoptosis) as follows: (1) group A, miRNAs with a significant up/downregulation at P13 and P15 compared with age-matched WT (P13 and P15) and rd10 at P13; and (2) group B, miRNAs with a differential expression at P17 compared with age-matched WT (P17) and with rd10 at P13: with no differential expression found at previous stages. In line with these criteria, we considered group A to contain the miRNAs that were most likely to be involved in the disease etiology and group B the miRNAs that were more likely to be involved in compensatory mechanisms and/or in disease progression. 

miRNA–mRNA Interaction Network Analysis and Gene Ontology and Pathway Enrichment Analysis

In silico prediction of mRNA targets for DE miRNAs was obtained using the miRWalk 2.0 database (available in the public domain, http://zmf.umm.uni-heidelberg.de/apps/zmf/mirwalk2/),³² which is able to predict targets for miRNAs using the algorithms of several databases including miRWalk 2.0 itself, miRanda,^33,34 TargetScan,³⁵ and RNA22,³⁶ among others. Using Cytoscape v3.0.0 open access software,³⁷ we constructed miRNA–mRNA interaction networks with both DE miRNAs and their target mRNAs as predicted by miRWalk 2.0. Only genes predicted by at least two algorithms were considered for further analysis. 

To identify inverse expression relationships among miRNAs and mRNAs, we filtered these networks by incorporating the data from our mouse transcriptome arrays, selecting mRNAs with enriched expression in the retina (signal values above background) and mRNAs with an opposite fold change to that of their predicted regulator miRNAs. We generated three interaction networks, namely group A miRNAs versus DE mRNAs at P15 and at P17 and group B miRNAs versus DE mRNAs at P17. 

Figure 1

 

Photoreceptor cell death in rd10 mouse retina. (A, B) PR cell death (green nuclei) starts abruptly at P16 in the central retina, with no apoptotic nuclei detected neither at P13.5 (data not shown) nor at P15.5. (C–G) Cell death progresses toward the periphery. (H) Schematic representation of the PR cell death process. Preapoptotic nuclei are represented by blue circles and apoptotic nuclei are represented by green circles. Red squares represent the retinal areas where images were captured. Cell nuclei (blue) are labeled with 4′,6-diamino-2-fenilindol (DAPI), GCL, ganglion cell layer; INL, inner nuclear layer; ONL, outer nuclear layer; ON, optic nerve.

Figure 2

 

miRNAs expression levels in rd10 mouse retina. Differential miRNA expression at postnatal days 13, 15, 17, 19, and 22 are measured as the fold-changes compared with age-matched WT samples (left) and rd10 P13 samples (right). (A, B) Group A miRNAs. (C, D) Group B miRNAs. See Materials and Methods for a description of groups A and B.

Figure 3

 

miRNA–mRNA interaction networks. A selection of representative differentially expressed miRNAs (red squares) and their predicted target genes with an inverse expression relationship (green circles). (A) miRNA–mRNA network of differentially expressed miRNAs prior to the onset of apoptosis (group A) and (B) just after the onset of apoptosis (group B). Arrows indicate direct interactions between an miRNA and its target genes. Changes in the shade of green of mRNA nodes reflect differences in fold-change, starting from 1.3-fold (light green) to up to −4.33-fold in group A and −7.65-fold in group B (dark green).

The differentially expressed putative target genes whose transcript levels were inversely expressed relative to the relevant DE miRNAs (inversely expressed genes) underwent analysis of gene ontology. For this purpose, we used CytoScape's ClueGO app (version 2.2.3) and analysis of biological pathways³⁸ based on the Database for Annotation, Visualization and Integrated Discovery (DAVID, available in the public domain, http://david.abcc.ncifcrf.gov/) and WikiPathways (available in the public domain, www.wikipathways.org) to identify biological processes affected by miRNA dysregulation. 

Validation of miRNA and mRNA Expression by Quantitative PCR

Expression levels of DE miRNAs and inversely expressed mRNAs were analyzed for validation using a second technique, quantitative PCR (RT-qPCR). Validation was performed in the same samples as those studied in array analysis (i.e., three samples for each rd10 and WT mouse at three different time points [P13, P15, and P17]). RT-qPCR was conducted using miScript Primer Assays and the miScript SYBR Green PCR Kit (Qiagen, Hilden, Germany) for miRNA expression and the SYBR Green PCR Master Mix (Applied Biosystems, Foster City, CA, USA) for mRNA expression, according to the manufacturer's protocol. Each biological sample was amplified in triplicate. Threshold cycle values for miRNA expression were normalized to the mean expression values of miR-26a and miR-191-5p, which showed lowest variability of expression levels among all samples analyzed. For normalizing mRNA expression, we used the mean values of Gapdh, Tubb5, and β-Act genes for mRNA expression, which showed strong correlation between themselves in their levels of expression across all samples (r² = 0.905 and P < 0.001). 

Results

Using the rd10 mouse model of RP, we conducted a detailed global study of the expression of more than 1900 miRNAs and the correlation of the levels of expression with the whole transcriptome to investigate possible differences in miRNA and mRNA expression profile between rd10 and WT mice retinas. This allowed us to identify a set of differentially expressed miRNAs and to analyze their possible role in the biological pathways involved in retinal degeneration through their predicted target genes. 

TUNEL Analysis

We carried out a detailed analysis of photoreceptor cell death dynamics in retinal sections obtained from mice ranging from P13 to P22 at 12-hour intervals to select the time points of interest for the miRNA–mRNA expression studies (Fig. 1). No TUNEL signal indicating apoptosis was detected prior to P16 (Figs. 1A–1F). A signal appeared on P16 in the central retina, mainly restricted to the outer nuclear layer (ONL) and progressed toward the peripheral retina, similar to what has been observed in previous studies.^28,29 Based on our results, we selected retinal samples at three key time points for miRNA and mRNA expression profiling studies: prior to the onset of apoptosis at P13 and P15 and just after the onset of apoptosis at P17. 

miRNA and mRNA Expression Profiling

Based on the photoreceptor cell death dynamics, we analyzed the overall expression of miRNAs and mRNAs immediately prior to and after the onset of apoptosis, at P15, and at P17, when virtually all photoreceptors were still present. We also analyzed rd10 retinas at P13 to profile miRNA expression at an early stage, 3 days prior the onset of apoptosis. We considered that any possible alteration in the miRNA profile that could be happening in the rd10 mouse at this early stage is less likely to be influenced, if at all, by the apoptotic processes at this time. 

The GeneChip miRNA 4.0 arrays revealed a large number of significantly differentially expressed (DE) miRNAs between the two strains with a greater than ±1.5-fold change (P < 0.05; Supplementary Table S1). Of more than 1900 miRNAs explored, 152 miRNAs were differentially expressed in the retina of rd10 mice with respect to control samples, and from this subset, we compiled a list of 19 candidate miRNAs for qRT-PCR validation analysis based on fold-change values and on preliminary prediction studies of their target mRNAs. In addition, miScript miRNA PCR arrays revealed a total of 18 DE miRNAs, this list being reduced to 6 candidate miRNAs by applying the same filtering criteria as for the GeneChip arrays (Supplementary Table S2). Therefore, our candidate miRNA list increased to 25 miRNAs, of which 19 were upregulated and 6 were downregulated. 

Validation of candidate miRNAs was performed by qRT-PCR. Of 25 miRNAs tested, 19 miRNAs were validated. Data were presented as fold change over controls (Table 1). Most of the validated miRNAs were upregulated (17 of 19 miRNAs), and only two of six downregulated miRNAs were validated by qRT-PCR. Among validated miRNAs, the most highly upregulated were miR-6240 and miR-6970 (22-fold and 21-fold at P17, respectively), whereas miR-20b-5p and miR-19b-3p were the only miRNAs that were downregulated (−1.76-fold and −1.73-fold at P17, respectively). 

Table 1

 

List of Differentially Expressed miRNAs Validated by RT-qPCR

 

The expression of three group B miRNAs, miR-155-5p, miR-142a-5p, and miR-146a-5p, were analyzed over a longer time frame, including measurements at P19 and P22. We selected these miRNAs based on the magnitude of the differential expression and/or the relevance of their predicted target genes in the disease. Expression of this group of miRNAs continued to increase at P19, reaching a plateau by P22, with expression values similar to those observed at P19 for miRs 155-5p and 146a-5p, and a marked drop in FC values from 12 to 3 for miR-142a-5p at P22 (Fig. 2). This pattern in expression levels from P19 might be indicating an exhaustion of compensatory mechanisms due to irreversible photoreceptor (PR) loss. Following the same selection criteria as that applied for group B, we selected four miRNAs to analyze over a longer time frame, namely, up to 22 days: miR-6240, miR-6937-5p, miR-3473b, and miR-7035-5p. The miRNA expression levels were still elevated at P19, decreasing to expression levels comparable to WT retinas by P22 (Fig. 2). 

Target Gene Analysis and miRNA–mRNA Networks

We identified an average of 900 target genes for each miRNA. Information from these target gene analyses was combined with data from miRNA expression profiling arrays to construct miRNA–mRNA regulatory networks using CytoScape open access software (Fig. 3; Supplementary Fig. S1). Then, we filtered these networks using data from our mouse transcriptome arrays. This approach allowed us to identify genes with levels of expression that were inversely correlated with those of our candidate miRNAs, which may play important roles in retinal degeneration. The 554 inversely expressed genes identified were analyzed using DAVID and WikiPathways open access databases to acquire knowledge about their biological relevance and to identify new actors and their involvement in the retinal degeneration process. 

Gene Ontology and Pathway Analysis

Target genes with an expression that was inversely correlated with that of the differentially expressed miRNAs subsequently underwent analysis of genetic ontology (GO) and biological pathways to improve our understanding of the mechanisms by which the candidate miRNAs might affect the rd10 mice phenotype. For this, we used DAVID and WikiPathways web services and CytoScape's ClueGO application. Ontology and biological pathway analysis revealed that a large number of inversely expressed targets are encoding components involved in biological functions likely to be involved in the pathogenesis of retinal degeneration, in compensatory mechanisms intended to prevent further retinal damage or in genes that are essential for retinal survival (Table 2). Considering these findings, we selected genes related to the pathways involved in apoptosis, inflammatory response, cytokine and chemokine activation, calcium signaling, visual perception, and visual phototransduction, among others, for qRT-PCR validation. In addition, we selected inversely expressed genes with high fold-change values, even if they had no clear involvement in retinal degeneration. 

Table 2

 

Summary of Differentially Expressed miRNAs and Their Predicted Target mRNAs

 

For target genes of the miRNAs that we classified into group A (i.e., differentially expressed prior to the onset of apoptosis), we mainly focused on pathways potentially involved in photoreceptor cell death such as apoptosis or inflammation, whereas for target genes of the miRNAs that we classified into the group B (i.e., differentially expressed just after the onset of apoptosis), we mainly focused on genes involved in pathways likely to be involved in photoreceptor survival. In addition to this, we also detected differentially expressed genes directly ascribed to retinal diseases such as Klhl23 or Nr2e3 involved in cone-rod dystrophies or RP, respectively. In total, we selected 49 genes for qRT-PCR validation. 

Expression levels of these 49 genes, 45 downregulated and 4 upregulated, were analyzed by qRT-PCR for validation. Downregulation of 36 (80%) and upregulation of 3 (75%) of the selected genes was confirmed, summing to a total of 42 genes that showed differences in expression between rd10 and WT retinal samples. Table 2 summarizes information about the fold regulation of the candidate miRNAs, their inversely expressed target genes and their associated signaling pathways that might be involved in the retinal disease. 

Finally, some of these genes have predicted interactions with more than one of the candidate miRNAs. For instance, Pde6g, which encodes the γ inhibitory subunit of the rod-specific phosphodiesterase 6, a key component of the phototransduction cascade, is a target for miR-3473b, miR-7035-5p, and miR-762. Another example is Peli3, a negative regulator of apoptosis, which is a predicted target for miR-6937-5p and miR-7115-5p. 

Discussion

In this work, we contributed to expanding our knowledge about the molecular mechanisms underlying RP. We identified a total of 19 candidate miRNAs that may play a fundamental role in the physiopathology and/or in the progression of the disease. GO and pathways enrichment analysis showed that some of the genes targeted by these miRNAs, which also showed differential expression between rd10 and WT mice, are involved in biological pathways with potential roles in retinal degeneration or survival. 

Given the established role played by miRNAs in retinal development and survival,^1,2 as well as in a growing number of distinct diseases,^{7–10,12,14,15,19,22,39,40} we hypothesized that alterations in miRNA expression might contribute to the molecular mechanisms involved in the physiopathology of RP. To provide experimental evidence to support the possible association between miRNA expression and retinal degeneration, miRNA expression profile in retinas from the rd10 RP mouse model^27–29 and C57BL/6J wild-type mice were compared by microarray and RT-qPCR analysis in the early stages of retinal degeneration. To gain further insight into the possible relationship between miRNA alteration and gene expression regulation with implications in retinal degeneration, we also performed transcriptome expression profiling of the same samples as used for miRNA expression analysis. 

We considered it likely that miRNAs down- or upregulated at a stage just prior to the onset of apoptosis might be associated with PR cell death mechanisms. Following the same rationale, we considered that miRNAs that are deregulated at latter stages might play relevant roles in compensatory mechanisms attempting to counteract PR death. Therefore, one of the key methodologic aspects of our study was to establish the time points for the miRNA–mRNA expression profiling analysis having first ascertained the exact time of onset of photoreceptor degeneration in our rd10 mouse model. 

According to previous studies, the onset of PR apoptosis in the rd10 mouse was reported to take place between the second and third postnatal weeks.^27–29,41 The onset of PR death is dependent on light conditions, with slight differences reported among different laboratories, depending on the intensity and duration of animals' exposure to light.⁴² To accurately establish the onset of the death of these cells, retinal samples were obtained every 12 hours between P13 and P22 and analyzed using TUNEL assay. We determined that apoptosis starts in the central retina between P15.5 and P16, with PR cell death then progressing toward the periphery (Figs. 1D–1F), in accordance with previous reports.^28,29,41 

One of the pathways that we focused on was inflammatory response; previous studies have evidenced that sustained chronic inflammatory process may contribute to the etiology of RP in rd10 mice and in other animal models of retinal degeneration.^43–46 Furthermore, these studies support the notion that miRNAs are rapidly upregulated in response to inflammatory signals and may either stimulate the magnitude and duration of inflammation or silence it.⁴⁷ Example of some of the genes related with inflammatory or immune response included Sema7a, Pan2, Pla2r1, Thbs1, or Trim36. Taking into account that deregulation calcium homeostasis plays a key role in retinal degeneration in rd10 mice, pathways involved the regulation, signaling, and/or exocytosis of this cation were also selected for validation, such as those related with Doc2b or Pde6g genes. In addition, we selected genes involved in physiologic processes or in the normal function of retinal cells, such as intraciliary transport (Lac5), phototransduction signal (Rgr), or retinol metabolism (Ttr). 

Interestingly, among the inversely expressed genes, there were genes previously linked to retinal disorders. For instance Klhl23 and Nr2e3, which are predicted targets of miR-6240, are linked with cone-rod dystrophy⁴⁸ or RP,⁴⁹ and Rrg, a predicted target of miR-3473b is associated with chorioretinal atrophy.⁵⁰ 

With respect to most of the putative target genes, they showed a modest differential expression, in comparison with that found in the expression of their associated miRNA. That is, despite some miRNAs showing fold-change values greater than 5 or even 10, their predicted target genes showed a significant, although rather modest changes in their expression levels of about 1.5- to 2-fold, with the exception of Serpina3n, Xaf1, and Osmr, with fold-change values of 113, 7, and 7, respectively (Table 2). This difference in the magnitude of the expression levels of miRNAs and their target genes could be a consequence of the complex regulatory network in which a relatively few miRNAs can regulate as many as 50% to 60% of genes in humans and other mammals. Therefore, convergence of miRNA action on hundreds of genes implies that single miRNAs may have limited effect on the expression of their mRNA targets.⁴ Moreover, it should also be considered that some of the candidate miRNAs might be inhibiting protein translation rather than degrading mRNA. In such cases, mRNA levels would not be affected. Notwithstanding the fact that they showed significant but relatively small differences in expression, they are of potential interest in the context of molecular mechanisms underlying RP. 

It is also remarkable that in most of genes with an inverse expression to that of group A miRNAs, that is, miRNAs deregulated at P15, showed significantly differential expression at P17, but not at P15. This delay might be related to the time that miRNAs require to exert their silencing effect on their mRNA targets. To the best of our knowledge, this is the first study that has described a time-delayed interaction between differentially expressed miRNAs and their predicted mRNA targets. This was possible through transcriptome analysis of retina samples at both P15 and P17, which allowed us to identify reverse expression relationships between DE miRNAs and their target mRNAs in a time resolved manner. In this regard, we considered that a 2-day lag time would allow a given miRNA to exert its inhibitory effect over its target mRNAs. Had we only studied miRNA and mRNA profiles at the same time points, the identification of this interaction would have been missed in most cases. 

There are previous reports of the involvement of miRNAs in retinal development, function, and photoreceptors survival^1,2 and also the involvement of different miRNAs in the pathogenesis of various diseases of the retina including miR-9, miR-34a, miR-125b, and miR-155 in macular degeneration,^17,18 miR-146a and miR-195 in diabetic retinopathy,^19,20 and miR-125a and miR-17 in retinoblastoma.^21,22 

Interestingly, some of the differentially expressed miRNAs we detected just after the onset of apoptosis are common to some of these retina-related miRNAs, including miR-146a, miR-155, miR-29b, and miR-17 (Table 1), as described below. 

MiR-17 is a member of the miR-17-92 cluster that has been implicated in inflammatory diseases,⁵¹ Burkitt lymphoma,⁵² and retinoblastoma.⁵³ Some evidence suggests that this cluster regulates a broad spectrum of biological processes of T-cell immunity including a negative expression relationship between miR-17 and the synthesis of Sirp-α gene that effectively regulate macrophage inflammatory responses.⁵⁴ 

Among the retina-related miRNAs described above, we only found miR-17 and miR-29b to be inversely expressed with genes involved in pathways with a putative involvement in retinal degeneration. Examples of this are Mtfp1 and Mllt11, both target genes for miR-17, and involved in positive regulation of apoptotic process.^55,56 With regard to miR-29b, we detected downregulation of Pla2r1, Pan2, and Duxbl2, genes targeted by this miR that have been implicated in inflammatory and immune responses.⁵⁷ Hence, increased expression of miR-17 and miR-29b could be part of survival mechanisms initiated to slow down the PR cell death process via the modulation of these genes. 

As far as RP is concerned, Loscher et al. reported upregulation of miR-1a, miR-133a, and miR-142a in four different mouse models of RP linked to genes involved in both autosomal dominant and autosomal recessive forms of the disease, rhodopsin and rds/peripherin, respectively. Our results in rd10 mice are consistent with their findings.^25,26 These results reinforce the hypothesis that there are commonalities in the miRNAs expression pattern in different types of RP regardless of the nature of their genetic cause, which could open the way to common therapeutic strategies applicable to different types of RP. 

However, the role that upregulation of either miR-1 alone or the miR-1/133a cluster might be playing in the apoptotic process is not exempt from controversy. There are data suggestive of both anti- and proapoptotic effects on myocardiocytes under ischemic conditions, by targeting of the proapoptotic genes Casp9 and Bcl2^58,59 and prosurvival genes including Pkce; Hsp60, and Igf1,^60,61 respectively. 

In any case, differential expression of miR-1a, miR-133a, and miR-142 in our study, as well as in the work of Loscher et al., is observed after the onset of apoptosis. This might be reflecting an activation of compensatory mechanisms in an attempt to prevent photoreceptor cell death in the mutant retina. In support of this, selective ablation of Müller cells resulted in photoreceptor death, and an altered expression of miR-1a, miR-133a, and miR-142.³⁹ Hence, the observed increase in the expression of miR-142a and the miR-1a/133a cluster might be a direct consequence of the cell death process and a common feature in diseases that cause inflammation and apoptosis. 

Among genes expressed inversely to miR-142a, we found genes involved in apoptosis (Trp53bp2 and Rraga), and in inflammatory responses (Thbs1 and Trim36). Trp53bp2 acts as a positive regulator of the execution phase of apoptosis by encoding a member of the ASPP family (apoptosis-stimulating protein of p53) and Rraga plays a direct role in a Tnfα signaling pathway leading to induction of cell death.⁶² On the other hand, Thbs1 is involved in inflammation through the activation of phospholipase C enzymes, which results in intracellular calcium levels rising and the activation of protein kinase C (PKC) in platelets,⁶³ and Trim36 is part of the superfamily of tripartite motif-containing (TRIM) proteins, which are induced by type I and II interferons and are crucial for induction of effective immunity.⁶⁴ The fact that the expression of these genes decreased in the early stages of retinal degeneration may be indicating an attempt to counteract the cell death process by increasing the expression of miR-142a. 

Differential expression of miR-7036b-3p, miR-7221-3p, and miR-762 was found at P17 compared with WT samples with no differential expression detected at P15. Among them, only one has previously been linked to disease, namely miR-762, which has been linked to vascular smooth muscle cell calcification⁶⁵ and breast cancer.⁶⁶ Notably, in the current study, we detected an inverse expression relationship between miR-762 and the Pde6g gene. As already discussed, Pde6g is a key player in the phototransduction cascade. Hence, increased levels of miR-762 might be contributing to exacerbate the death of rod cells in this mouse model with mutations in the Pde6b subunit by silencing the expression of the inhibitory γ subunit Pde6g of the PDE6 protein complex. 

In addition, we found an inverse expression relationship of miR-762 with Rgr, Trafd1, Nnat, E2f3, and Otud7b mRNA, whose downregulation might be involved in the physiopathology of retinal degeneration, these genes having roles in relevant retinal processes such as phototransduction, negative regulation of the innate immune response, maintenance of nervous system structure, and regulation of transcription and immune system processes, respectively, according to GO and pathway enrichment analysis. 

Interestingly, we found decreased expression of putative target genes for miR-7036b-3p and miR-7221-3p that might be involved in the development of retinal dystrophy rather than in compensatory mechanisms. Specifically, the target genes for miR-7036b include Ptn, a negative regulator of neuron apoptosis, and Crx, a causal gene of RP and other inherited retinal dystrophies. Further, one predicted target gene for 7221-3p is Rdh8, which encodes a retinol dehydrogenase that catalyzes the reduction of all-trans-retinal to all-trans-retinol, the first reaction step of the rhodopsin regeneration pathway, and its deficiency has been linked to AMD.⁶⁷ 

Nevertheless, although it is possible that the above-mentioned genes, with significantly reduced levels of mRNA expression, are involved in the retinal degeneration in the rd10 mice, differential downregulation of miR-762, miR-7036b-3p, and miR-7221-3p takes place after the onset of PR apoptosis. Therefore, further studies are required to discern whether these three miRNAs are contributing to the physiopathology or to the progression of the retinal dystrophy. 

Group B miRNAs, which have been previously associated with retinal diseases, were up- or downregulated after the onset of retinal degeneration. This might be indicating an involvement of these miRNAs in compensatory mechanisms against PR cells death or a deregulation resulting from PR cell death itself, and it is considered less likely that they contribute to disease onset. 

In contrast, group A miRNAs showed differential expression at P15, prior to the onset of apoptosis. Therefore, this group of miRNAs may play important roles in the etiology of retinal dystrophy. Strikingly, to the best of our knowledge, none of the miRNAs included in group A have been previously linked to any retinal disease. 

Group A is composed of seven miRNAs: miR-6937-5p, miR-6240, miR-3473b, miR-6970-5p, miR-7035-5p, miR-7081-5p, and miR-7115-5p. Notably, overexpression of all these miRNAs had increased by P17. For instance, miR-6937-5p showed a 2-fold upregulation at P15 and a 10-fold increase at P17. 

Because most of the miRNAs we include in group A have only been described in recent years, there is little information regarding their possible involvement in different biological pathways. To the best of our knowledge, these miRNAs have not been ascribed to any disease. The expression of miR-6240-5p was more than 22-fold higher in rd10 retinas, this being the miRNA that showed the greatest upregulation. Among the target genes for miR-6240-5p, we found two genes associated with RP, Klhl23 and Nr2e3, and one with retinoschisis, Rs1. Another noteworthy target for miR-6240 was Nxnl2, which has been suggested to be involved in the maintenance of both the function and the viability of sensory neurons, including photoreceptors.^68,69 

The second most upregulated miRNA was miR-6970-5p (Table 1). Inverse expression studies revealed Eif2s3y, Onecut1, and Tspan6 genes to be negatively expressed with this miRNA. Eif2s3y participates in the GTP catabolic process⁷⁰; Onecut1 in transcription and cell cycle regulation⁷¹; and Tspan6 mediates immune signaling in a ubiquitination-dependent manner.⁷² 

miR-6937-5p showed a significant 10-fold increase at P17. Its predicted targets, Doc2b, Lca5, and Peli3, were observed to be downregulated in our experiments. Doc2b is expressed ubiquitously and is suggested to be involved in calcium-dependent intracellular vesicle trafficking in various types of cells.^73,74 Diseases associated with the Lca5 gene include Leber congenital amaurosis and severe early-childhood-onset retinal dystrophy.^75,76 This gene has functions in photoreceptor cell maintenance and in intraciliary protein transport as part of the centriole.⁷⁷ With regard to Peli3, this gene acts as a negative regulator of the extrinsic apoptotic signaling pathway, and its suppression leads to enhanced formation of the death-induced signaling complex in response to TNF.⁷⁸ 

Two interesting putative target genes for miR-3473b that we found to be downregulated were Rgr and Pde6g,^79–81 both involved in the phototransduction signaling cascade: Rgr as part of the G-protein–coupled receptor signaling pathway and Pde6g encoding the gamma subunit of cyclic GMP-phosphodiesterase. Pde6g expression is restricted to PR, whereas Rgr is expressed in PR, RPE, and Müller cells. 

Taking into account all the aforementioned findings regarding group A miRNAs, including the observed differential expression prior to the onset of apoptosis, the relevant role of their target genes for the retina in health and disease, and the sustained increased expression until at least P19, our data support the view that it is likely that this group of miRNAs plays an important role in the etiology of retinal dystrophy. 

The observed downregulation/suppression of genes involved in apoptosis inhibition and inflammation regulation is in good agreement with the apoptotic death of PR cells and inflammatory processes observed in rd10 retinas, emphasizing the role that differentially expressed prior to the onset of apoptosis miRNAs might be playing in the retinal etiology. 

Although it was initially proposed that miRNAs act mainly at the protein translation level, it is now widely accepted that both mRNA degradation and translation repression occur. However, despite great progress in elucidating the mechanism of miRNA-mediated mRNA targeting has been made, the extent and frequency of these two mechanisms has not yet been determined.^82,83 In this regard, one of the limitations of our study is that it is focused only on the analysis of mRNA expression level. 

Another limitation has to do with the predictive nature of our study. Despite many sequence-based target computational prediction algorithms having been developed and widely used to predict miRNA targets, their use implies inherent error in target predictions.⁸⁴ These computational tools rely not only on the sequence complementary, but on thermodynamics, energy potentials, accessibility of target sites, and conservation, essentially scanning the entire transcriptome to predict a set of targets for each miRNA.⁸⁵ However, this has the limitation of not considering that the relationship between miRNA and target genes is dynamic. This means that under different conditions, the same miRNA could have different sets of targets, and therefore different functions, when expressed in different cell types. Thus, these methods produce many false positives, even if the target information is accurate.⁸⁴ To control false positives, a strict cutoff value is typically established; however, in this case, these methods tend to reject many true target relationships, that is, false negatives. Therefore, in this target prediction study, it is necessary to bear in mind that there is an inherent error that could condition our posterior results and that the inverse expression of the target mRNA of a particular miRNA is not at all proof of miRNA–mRNA regulation. Thus, miRNA–mRNA target pairs should be further validated in biological experiments, and therefore what our study is providing are strong candidates, which need further experimental validation. 

Conclusions

In the present work, we identified 19 DE miRNAs and 39 inversely expressed target mRNAs that may be key players during retinal degeneration, considering the reported involvement of these DE miRNAs and mRNAs in basic retinal processes such as apoptosis and inflammation regulation, which take place in degenerating rd10 retinas. 

Among the miRNAs identified, seven have previously been reported in retinal studies, whereas nine have not previously been ascribed to any diseases. 

As far as we are aware of, this is the first study describing a time-delayed interaction between differentially expressed miRNAs and their predicted mRNA targets. The identification of disease-related miRNAs and the consequent alteration in gene expression will contribute toward our understanding of the etiology and evolution of retinal degeneration. 

Acknowledgments

The authors thank Carlos San José, Elisabeth Hijona, Mariasun Zabala, Iñaki Martinez, and Manuel Yecora at the animal facility at Biodonostia Research Institute for help in the maintenance of the mice colony. 

Supported by the Fundación Jesús de Gangoiti Barrera and from the Basque Government's Department of Industry and Education Grants SAIOTEK-PE11BN002, PC12BN001, and DEPLC13/002 (AA, JRE); funds from Foundation of Patients of Retinitis Pigmentosa of Gipuzkoa (Retinosis Gipuzkoa Begisare) (OB); a grant from the Fundación Mutua Madrileña (OB); Basque Government's Department of Education grants DEDUC14/309 (MEI), Diputación Foral de Gipuzkoa DFG15/006 (MM-C), and ELKARTEK 16/014 (MM-C); National Institute of Health Carlos III (Instituto de Salud Carlos III) Grants ISCIII: CP10/00572 (JRE), PI13/02621 (JRE); an Intensificación Contract (ALdM) from the Basque Government's Department of Industry; and a grant from the Foundation of Patients of Retinitis Pigmentosa of Gipuzkoa (Retinosis Gipuzkoa Begisare) (JRE). JR-E is a Miguel Servet II Fellow, National Institute of Health Carlos III (Instituto de Salud Carlos III), ISCIII: CPII16/00012. 

Disclosure: A. Anasagasti, None; M. Ezquerra-Inchausti, None; O. Barandika, None; M. Muñoz-Culla, None; M.M. Caffarel, None; D. Otaegui, None; A. López de Munain, None; J. Ruiz-Ederra, None